MOLECULAR BASED DECISION SUPPORT SYSTEM FOR CANCER

TREATMENT

TECHNOLOGY FIELD The present invention relates to a method of generating treatment prioritized recommendations for cancer patients, based on molecular data and statistical Evidence-Based Medicine.

CROSS-REFERENCE TO RELATED PATENT APPLICATIONS This patent application claims priority from and is related to U.S. Provisional Patent Application Serial Number 62/007,938, filed 5 June 2014, this U.S. Provisional Patent Application incorporated by reference in its entirety herein.

BACKGROUND Molecular Medicine strives to understand normal body functioning and disease pathogenesis at the molecular level, which may allow researchers and physician-scientists to use that knowledge in the design of specific molecular tools for disease diagnosis, treatment, prognosis, and prevention.

A considerable amount of molecular data relating to a cancerous tissue is available to physicians today, such as:

- HIC panel / Expression profiling

- Full chemo-sensitivity assay

- "Hot-spots" RT-PCR panel

- "Deep" / NG sequencing

- "Liquid-biopsy" study (CTCs / n/t-DNA) Currently, oncologists have three main sources of data to help them recommend a specific treatment to a cancer patient:

- The patient's medical records, including cancer type, stage of disease, previous treatments, etc.

- Molecular data derived from any of the above assays, which may

provide a table, such as depicted in Fig.1 , comprising, for each marker tested:

• The associated treatment

• The molecular assay type, strength and extent

· An indicator showing whether the patient is positive (has a

potential to react to the associated treatment) or negative.

- Publications of clinical researches showing statistical distributions of the various treatments efficiency among sampled populations. The efficiency is measured in terms of overall survival (OS), remission time (static disease included) (TTP) and response rate (RR), such as depicted IN Fig. 2.

The main drawback of fitting a treatment to a patient using only these available data sources lies in the fact that there is no known way of comparing / prioritizing a number of positively indicated treatments. The clinical researches are done on molecularly heterogeneous populations, which are sampled according to various parameters which do not include molecular homogeneity, but rather include other parameters such as age, gender, cancer type, tumor stage etc., which are less relevant when it comes to cancerous tumors. Even if molecular homogeneity is included in the sampling parameters - the tumor heterogeneity is by far more robust for a clinical trial to contain.

There is need for a method that will predict the probability of a patient to benefit from a treatment and the benefit itself, given that the patient has a positive marker for this treatment. SUMMARY

According to an aspect of the present invention there is provided a method of determining the statistical distribution of a treatment's benefit to a patient having a positive molecular marker to said treatment, comprising: using general population statistics to determine the positive effect of said treatment on a patient; using general population statistics to determine the conditional probability of a patient with positive marker to react to said treatment; and using said determined positive effect and conditional probability to determine a benefit factor of said treatment to the patient. Determining the positive effect may comprise estimating the survival function of patients reacting to the treatment using general population survival function.

According to another aspect of the present invention there is provided a molecular based decision support method for cancer treatment, comprising: performing the method of claim 1 for a plurality of treatments for which a patient has positive molecular markers; and ranking said treatments according to said determined benefit factors.

Ranking may comprise calculating the significance of the difference between efficiencies of the various treatments and determining a resolution cutoff rule for differentiating between various treatments. According to yet another aspect of the present invention there is provided a system for determining the statistical distribution of a treatment's benefit to a patient having a positive molecular marker to said treatment, comprising: a system server running molecular guided analysis; at least one patients history repository connected with said system server; at least one medical publications repository connected with said system server; and at least one source of molecular laboratory test results connected with said system server, said system configured to supply personally prioritized treatment recommendations to a patient. BRIEF DESCRIPTION OF THE DRAWINGS

For better understanding of the invention and to show how the same may be carried into effect, reference will now be made, purely by way of example, to the accompanying drawings.

With specific reference now to the drawings in detail, it is stressed that the particulars shown are by way of example and for purposes of illustrative discussion of the preferred embodiments of the present invention only, and are presented in the cause of providing what is believed to be the most useful and readily understood description of the principles and conceptual aspects of the invention. In this regard, no attempt is made to show structural details of the invention in more detail than is necessary for a fundamental understanding of the invention, the description taken with the drawings making apparent to those skilled in the art how the several forms of the invention may be embodied in practice. In the accompanying drawings:

Fig. 1 is a table showing molecular data derived from various assays; Fig. 2 is a table showing efficiency measurements of various treatments;

Fig. 3 is a schematic drawing showing the various components participating in the system of the present invention; Fig. 4 is a general flowchart showing the main processes performed by the molecular guided analysis system;

Fig. 5 shows division of the Kaplan-Meier curve calculated for a given population;

Fig. 6 is a graphic representation of the overall survival; Fig. 7 shows estimation of RTD from Kaplan-Meier curve; and

Fig. 8 represents the step function for determining the number of different treatments. DETAILED DESCRIPTION OF THE INVENTION

In the following detailed description, numerous specific details are set forth regarding the method and the environment in which the method may operate, etc., in order to provide a thorough understanding of the present invention. It will be apparent, however, to one skilled in the art that the present invention may be practiced without such specific details. In other instances, well-known components, structures and techniques have not been shown in detail to avoid unnecessarily obscuring the subject matter of the present invention. Moreover, various examples are provided to explain the operation of the present invention. It should be understood that these examples are exemplary. It is contemplated that there are other methods and systems that are within the scope of the present invention. Also, the same reference numerals are used in the drawings and in the description to refer to the same elements to simplify the description. The present invention attempts to overcome the shortcomings of existing treatment recommendation methods by deducing, from existing PAB (Population Average Based) statistics, the statistical distribution of a treatment's benefit to a patient having a positive molecular marker to this treatment. Fig. 3 is a schematic drawing showing the various components participating in the system 300 of the present invention, comprising at least one patients history repository 310, at least one repository of medical publications 320 downloaded from the internet 325 or retrieved from local databases 327 and at least one source of molecular laboratory test results 330, all contributing to the molecular guided analysis system 340 for supplying personally prioritized treatment recommendations 350 to a patient.

Fig. 4 is a general flowchart 400 showing the main processes performed by the molecular guided analysis system 340.

In step 410 the system determines to which treatments the patient has a potential to react (positive marker) according to the molecular assays results. In step 420 medical publications relating to these treatments are retrieved from local databases 327 or by ad-hoc internet 325 searches. In step 430 general statistics is used to determine the positive effect of each treatment on the patient, as will be explained in detail below. In step 440 general statistics is used to determine probability of a patient with positive marker to react to each one of the treatment, as will be explained in detail below. In step 450 the general benefit of each treatment is calculated using the previously calculated values and in step 460 the various treatments are prioritized according to the calculated general benefit. Determining the positive effect of a treatment is done using the Kaplan-Meier estimator. The Kaplan-Meier estimator, also known as the product limit estimator, is an estimator for estimating the survival function from lifetime data. In medical research, it is often used to measure the fraction of patients living for a certain amount of time after treatment. A plot of the Kaplan-Meier estimate of the survival function is a series of horizontal steps of declining magnitude which, when a large enough sample is taken, approaches the true survival function for that population. The value of the survival function between successive distinct sampled observations ("clicks") is assumed to be constant. As shown in Fig. 5, the method of the present invention divides the Kaplan- Meier curve calculated for a given population, having an exponent fitting curve 510 and a hazard function lambda (e.g. 32.5), into two curves:

- A first curve having an exponent fitting curve 520 and a hazard function lambda (e.g. 10), showing survival rate for responders to the treatment (positive marker).

- A second curve having an exponent fitting curve 530 and a hazard function lambda (e.g. 50), showing survival rate for non-responders to the treatment (negative marker).

The hazard function lambda is defined as the event rate at time t conditional on survival until time t or later (that is, T > t). Survival function calculation

In the following text, parameter's notation in cursive capital refers to "true parameters", and in standard capital to their estimate; e.g. OS is the median Overall Survival of the study's sample and it's the estimate of an unknown "true" median that we will note os .

There is a probability p that the patient will respond after the treatment and a probability (1 -p) that the patient will not respond.

Let's call the survival function for patients who respond f0(x) and f1 (x) the survival function for patients who don't respond. Thus, the overall survival function is the mixture:

S(x) = p. f0(x) + (1 -p). f1 (x) (1.1 )

(f1 (x) may be different from survival function of the disease without a treatment because of side effects of the treatment)

We assumed that this relation is true for the medians (It's indeed true for the means but means are less usual in survival analysis, and simulations show that the bias for the mixture medians is not so substantial):

OS = p.RS + (1 -p).NRS (1 .2)

OS (Overall Survival) is the median time from the treatment to death/loss of follow up. RS (Response Survival) is the median time from response to death/loss of follow up.

NRS (Non Response Survival) is the median time from "non-response" to death/loss of follow up. In oncology p is called RR (Response Rate) and represents patients whose cancer shrinks (termed a partial response, PR) or disappears after the treatment (termed a complete response, CR).

In studies, OS and RR are frequently reported but not RS and NRS. However there is another parameter which is frequently reported in oncology studies and which may help reproduce RS and NRS:

The TTP (Time To Progression) or PFS (Progression Free Survival) is the median time from the treatment to the progression of the disease.

If we change the definition of "response" in the parameters RR, RS and NRS so that they also include patients with stable disease after the treatment (SD) (and only these ones who have Progressive Disease (PD) are excluded): we can assume (See Fig. 6) :

RS = TTP + NRS (1 .3)

OS = RRxRS + (1 -RR)xNRS (1 .4)

Hence:

(1 .5)

= ΊΤ<Ρ + OS - ⁽ R x ΤΤΡ = OS + (l - ¾ )x ^ (1 .6)

Where , os and Τ<ρ are provided from literature.

RR from the molecular results and its distribution

We now attempt to calculate the conditional distribution of RR, namely the probability of response given an existing marker.

RR is generally given in the study and is a proportion, so n*RR is binomial distributed: n.RR ~ Binom (RR , n) and we get:

RR.(l - RR)

Var{RR) (2.1 )

n

So we can compute

RR.(l - RR)

sd(RR)= (2.2)

n

RR.(l - RR)

So asymptotically, we get RR ~ N( RR, ) (2.3)

n

When we don't know about molecular markers significant results on effectiveness of a treatment (Molecular Prediction of Response = 0) we'll use (2.6) for estimating RR.

If we know that a treatment cannot work (Molecular Prediction of Response -1 ), we'll admit in a deterministic way that

If we know that a present marker promotes the treatment we can distinguish three cases:

First, let's define the events "the marker exists" Λ+ ,"the marker doesn't exist" Λ- and the events "the treatment works (including stable disease) " Ω+ ,"the treatment doesn't work" Ω-.

We make the biologically defendable assumption

Λ- > Ω- (or Ω+ > Λ+) (2.5)

Namely, if the marker doesn't exist the treatment doesn't work (or if the treatment works the marker exists).

• If the article provides the conditional probability of the treatment

working given that the marker exists Ρ(Ω+|Λ+) we'll use it and its sd is calculable in the same way as normal RR (See equation 2.2) • If we have relevant knowledge on the probability of marker existing Ρ(Λ+) we ΊΙ use a Bayesian estimation

Ρ(Ω+|Λ+) = RR/ Ρ(Λ+) (2.6)

• If we have neither, we'll use in a deterministic way

Ρ(Ω+|Λ+) = RR + a.sd(RR) (2.7) with a=3 for chemotherapy treatments and 2 for other treatments.

The table below summarizes the estimation of RR for the various cases:

Summary for estimating the Molecular Prediction of Response (MPR) parameter RR :

1 0

The article Ρ(Ω+|Λ+) Ρ(Ω+|Λ-) gives

sd calculable sd

Ρ(Ω+|Λ+) calculable other articles RR/ Ρ(Λ+)

Knowledge

give

from Without sd RR , sd(RR)

articles Ρ(Λ+)

sd(RR) no knowledge RR+a.sd(RR) ae{2,3} Without sd

The distribution of TTP, RS, NRS and ART, ARS, APS

We are interested in comparing TTP and RS and NRS from different studies (each study is a different treatment which may concern the patient), but these parameters only make sense when they are adjusted by their probabilities of benefit to the patient. Let's define new parameters ART, ARS and APS which signify the probability distribution of the various observed values of TTP, RS and NRS:

ART = RRxTTP (3.1 ) The probability of a given TTP value given positive reaction to the treatment. ARS = RRxRS (3.2)

The probability of a given RS value given positive reaction to the treatment.

APS = (1 -RR)xNRS (3.3)

The probability of a given NRS value given negative reaction to the treatment.

We assumed that the estimates of TTP, RS and NRS are normally distributed and we can estimate their standard deviations from a confidences interval that appears in the study.

But in survival analysis, software (excepting SPSS) and studies generally use the method of Brookmeyer and Crowley to compute a Confidence Interval (CI) for median overall survival. And this method doesn't use standard deviation at all and gives an asymmetrical CI.

But, as a median, we have:

where f is the unknown TTP's PDF (probability distribution function) and n is its sample size. To be conservative, we would choose the biggest one half CI (generally the upper one), to estimate the standard deviation, but the methods used by studies, previously conservative (discrete tests dual to CI)

So we defined (by calling henceforth TTP its estimate) upCI^TTP)- IOWCI(TTP)

sd(TTP)-- (3.5)

2 z

Where Z is a standard score.

In the case CI is not mentioned in the paper, we assumed that TTP fits an exponential distribution with median = Ιη(0.5)/λ so f(TTP)=0.5 and

asymptotically:

For RS and NRS we have their estimates from (1 .5) and (1 .6) and their standard deviations:

Var(NRS)= Var(pS - ART)= Var(pS)+ Var(ART)^> sd(NRS)=

-dVar(OS)+ Var{ART)

Var[( - RR)x TTP]= Var[RR x TTP]= Var(ARl) (3.7)

Var(RS)= Var(pS + (l - RR)x TTP)= Var(pS)+ Var(f - RR)x TTP)= Var(pS)+ Var(ART)= Var(NRS)

To find the standard deviation of ART = RRxTTP

we, conservatively, assumed that RR and TTP are independent and so: (3.9) And we get a table with the following format:

Articles/Treatments TTP Sd(TTP) RR Sd(RR) ART Sd(ART)

Article 1

We can create similar tables for ARS and APS, but in these cases, we don't use twice the sd of RR so we use it, the second time, in a

deterministic way: e.g. sd(ARS) = RR ² JVarpS)+Var(ART) (3, 10)

Treatment Prioritization

There are k possible treatments for the patient and for each treatment we have beforehand estimated by meta-analysis three ^AEV (Adjusted

Efficiency Values), namely:

ART = MPR x RTD, effOS = MPR x RS + (\ - MPR) x NRS, subjOS = ARS or APS (The description of these ^AEV and the model is found in appendix 1 below). Similarly, we have already estimated the standard errors of the different estimators assuming a normal distribution and conservatively assuming independence across them.

The aim of this algorithm is not to provide objective rating of treatments but only to rank them using relative grades: The highest score for the most efficient treatment will always be 10 and the lowest score will always be 1 (Naturally, this doesn't allow comparing two treatments of two different reports).

The "Total Efficiency Score" is a weighted average of 3 ^AEV grades with the following weights: ART 50%

effOS 45%

subjOS 5%

In this main text we'll describe how the algorithm rates the 3AEV . Since the procedure is exactly the same for each, we'll describe it once using the generic notation AEV , while the program executes it three times for ART , ARS and AP S .

Trivially, we could maintain the original ranking of the AEV estimators regardless their standard errors, i.e. a-priory to assume that the efficacy of each treatment is significantly different from the efficacy of any other. The confidence level at which the whole ranking is true is bound below by 1 /(k!) (probability of random guessing). This bound is very small and we can still improve it (see below), but in our work, bounding accurate probabilities by considering statistical hypothesis testing for the difference among any two treatment efficiencies, would be a very challenging task to perform. So we have chosen to consider only the vector of "Confidence Level of the Ranking" for each pair of successive treatment efficiencies (denoted from now on as "treatments", bearing in mind that we refer to their efficiencies scores). Let's denote the vector of probabilities of having efficiencies difference (which reflects the probability of measuring the efficiency of treatment "A" to be larger than the efficiency of treatment "B", given the averaged efficiencies and the standard deviations of both treatments) among any adjacent treatments by "

CLR " Looking at its basic properties, it is clear that each of the vector's components is bound below by 50%, as a result of random guessing (we will see later, how to compute it and we shall use this primary ranking to initialize the algorithm only). On this "random guessing" limit, the low confidence level allows the ranking to be well defined (however an intolerable confidence of this ranking is involved), reflecting k, well ranked, treatment groups of order 1 . The prediction regarding the rank of the treatments efficiency, however well defined, is useless because its lack of any statistical significance.

Looking at the opposite limit, it could be argued that any two treatments are absolutely different (statistically-wise), if and only if their CLR is 1 , hence 1 is an upper limit of those vector components. Imposing a very high confidence demand on a set of treatments, results in lack of any statistical valid difference between those treatments, or, in other words, having one (unranked) treatment group of order k. Such a predication is also useless to anyone. Between those extreme limits of CLR, we would want to split the treatments to some different "treatment groups" of one or more treatments in each, respecting the original ranking while also respecting the significance of the standard errors, so that we'll get a vector of Y treatment groups (2 < Y < k-1 ) and within each treatment group the difference between the different treatments is not statistically significant (in the sense that will be presented here below).

The purpose of the algorithm is to find a "not-too-small" number of treatment groups, which allows a "not too small" average CLR . If high resolution is required for the ranking (i.e. a large number of treatments groups), the average CLR will be decreased. Similarly, the number of treatment groups will be decreased if high reliability is required for the ranking (i.e. a large average CLR ). We use the term "resolution" referring to the number of treatment groups and we use the term "reliability" referring to the averageCZJ? . The value of the resolution parameter is defined over the set ,3,...,k] since a single group of treatments is not a pertinent result. We consider the

"resolution's range" by the interval [2,£], and similarly we define below a relevant "reliability's range". Given a set of treatment efficiency values, the aim of our algorithm is to find a CLR range such that both the resolution's range and the reliability's range are maximal. Initializing the algorithm (or step (0)):

The algorithm, as a whole, is illustrated by an example in Appendix II. We have an ordered vector of k AEV estimators

(AEV^ AEV^AEV^)

And its standard-errors vector: (We have also the corresponding sample sizes vector nf ^>' )

At each step, we also calculate AEG (Adjusted Efficiency Grades) with its standard-error :

AEG is a vector of grades such that to the higher AEV component the equivalent AEG component is 10 and to the lower AEV component the equivalent AEG component is 1 .

9 - {AEV _{I) - AEV^

AEG, + 1 (1.1)

⁽ 0 AEV _{K) - AEV _{{L) ( Λ} 9 - ά{ΑΕν ₍ Λ

6(AEG ₍ = ^ ⁰⁾⁷ 1.2

AEV _{K) - AEV _{1) where (1 .1 ) is the projection on the interval [1 , 10] and (1 .2) comes from

For each pair of successive (^EV^, AEV^%) we define δ$ = AEV % - AEV and we compute the probability to find a difference between the AEV estimators larger than or equal to <5^when the true difference is null. This probability is naturally the lower bound of the probability to find a difference larger than or equal to^as a function of the true difference when this true difference is negative or null (i.e. when the ranking is false). In classical hypothesis testing, for the null hypothesis "the ranking is false" this probability is called the confidence level of the test. Those probabilities form a vector of order k-1 which we denote as "Confidence Level of the Ranking": CLR ^

Assuming normal distribution of the estimators, with standard deviation parameters known to be exactly the standard errors estimated, we have:

CLR ≡ Φ > 0.5

Where≤Φ denotes the CDF (Cumulative Distribution Function) of the standard normal distribution:

((WWee ddoonn''tt uussee tthhee tt ddiissttrriibbuuttiioonn bbeeccaauussee ooff tthhee wwiiddee vvaarriieettyy ooff ssttaannddaarrdd eerrrroorrss 55 wwhhiicchh ccoommeess ffrroomm tthhee wwiiddee vvaarriieettyy ooff ssaammppllee ssiizzeess aammoonngg tthhee ddiiffffeerreenntt

ttrreeaattmmeennttss.. TThhiiss zz--tteesstt ffoorr ccoommppaarriissoonn ooff mmeeaannss iiss sspprreeaadd wwhheenn tthheerree aarree hheetteerroosscceeddaassttiicciittyy aanndd llaarrggee ssaammppllee ssiizzee))

TThheenn wwee rreettaaiinn tthhrreeee vvaalluueess::

•• tthhee ssmmaalllleesstt ccoommppoonneenntt ooff CCLLRR ⁽⁽⁰⁰⁾⁾ ccaalllleedd ""ccuuttooffff CCLLRR ⁽⁽⁰⁰⁾⁾ "" aanndd ddeennootteedd 1100 CCLLRR$;;

•• tthhee iinnddeexx ooff tthhiiss ccoommppoonneenntt ££ ⁽⁽⁰⁰⁾⁾ :: CCLLRR$$ == CCLLRR ;;

•• tthhee mmeeaann ooff tthhee CCLLRR^^''ss ccoommppoonneennttss,, ccaalllleedd ""aavveerraaggee CCLLRR ^^"" aanndd ddeennootteedd ξξ ^{{00}} ;;

CCLLRR yyiieellddss aa tthhrreesshhoolldd ffoorr tthhee rreelliiaabbiilliittyy ooff tthhee wwhhoollee rraannkkiinngg::

1155 ⁽⁽⁰⁰⁾⁾ IIss tthhee ccoonnffiiddeennccee lleevveell ooff ""ttrruuee rraannkkiinngg"" ffoorr aa ppaaiirr ooff ssuucccceessssiivvee AAEEVV eessttiimmaattoorrss rraannddoommllyy sseelleecctteedd..

TThhee oouuttppuutt ooff tthhee iinniittiiaall sstteepp iiss tthhaatt:: ""FFoorr aaCCLLRR ''ss tthhrreesshhoolldd iinn tthhee iinntteerrvvaall

treatment groups is hand the average CZR is ξ ⁽ ' ." 20 The general step (i) of the algorithm

At each step j , l < j < k - 2 , we represent two successive "old" treatments (or treatment groups) which bear some different efficiency values in the j-1 step, in one "new" treatment efficiency vales in the j step. This representation is described below (we refer to it by noting that we "pull" treatments' efficiencies 25 and represent it as a different efficiency value): • We are pooling (AEV^yAEV^^ as if both treatments were a single

From e sizes ⁿ ),

The first expression is a simple weighted mean, and the second one comes from the famous decomposition of the empirical variance: σ ² ( )= ^ ' - x ² n true for any sample of a random variable X .

That yields for two samples (x _l .., x _n ), ( ^...,^ ) for which only the descriptive empirical statistics are known: x,y,a _x ,a _Y \ e possibility to restore the standard deviation of the samples' union by:

∑ , ^{2 +} ∑^, ² fnx + myY _ n(p _x ² + x ² ) m^ +J? ² ) fnx + myY

n + m n + m J n + m n + m J Obviously we obtain the new sample size: n^ = + n^ .

The results of this pooling are three new vectors of k - j components. These "new" vectors are of order which is smaller by 1 compared to the "old" vectors. • We compute the current step CLR by the method mentioned above,

• We retain the three following values : the cutoff CLRfo' , its index ^J) and the average CLR ^(j) ^ ^iJ) We conclude for the step j ·.

"For aCLR 's threshold in the interval the (maximal) number of treatments groups is k - j and the average CZR is ξ ^Κ"

Remark about the pooling:

We decided to consider only the successive pairs of treatments and not all the pairs because of the algorithm's pooling: Indeed, it may happen that for three successive treatments groups A, B and C the confidence levels of the ranking between A and B and between B and C are larger than the confidence level of the ranking between A and C. In such a case, there is no interest to pool A and C without pooling them with B. Conclusion and final output

The algorithm described two steps functions of the cutoff's intervals. One is the average CLR function overall increasing and the second is the number of treatments groups which is decreasing. The intersection of both functions defines the optimal step J of the algorithm. It may be computed as already stated above, like the maximum of the product function. The product function is the steps function "product of the values on the same scale of both functions".

To project the two functions on the same scale we use the ranges [2,£]for the number of treatments and the range ξ ^{0 max ξ ^ω for the average CLR

1< j≤k-2

function. When the optimal step is not the same for the three AEG we calculate the overall optimum by the weighted mean of the 3 J found. i.e. J = round (0.5 x J^ + 0.45 x J _effOS + 0.05 x J _subjOS )

We retain also initial AEG and final AEG (the last step which is "necessary"). But to finally give the Total Efficiency Score ( TES ) we use Vne AEG ^(J) : which is a vector for k - J treatments groups.

Appendix I

Usually, in oncological studies we distinguish three reactions to a treatment well defined as their proportions are called:

ORR (Objective Response Rate which includes Partial Response and

Complete Response) based on a percentage reduction of the tumor, PD (Progressive Disease) based on a percentage increase of the tumor and SD (Stable Disease) in which there is neither sufficient shrinkage to qualify for response, nor sufficient increase to qualify for progression.

Our analysis is based on a two group mixture and we called the "responders group" the participants whose tumors did not progress. Their proportion equals to ORR+SD and is called CBR (Clinical Benefit Rate).

From a clinical trial or meta-analysis of several clinical trials on samples, which requires adequate similarities to our patient, we can estimate a naive probability p for our patient to respond (hereafter, the words "respond" and "response" are based on our definition and include stable disease) to the treatment (and a probability (1 -p) that the patient will not respond). The estimation is naive in the sense that we don't use any molecular knowledge specific to our patient. p = CBR Let's call the survival function for patients who respond f ₀ (x) and the survival function for patients who don't respond f-i(x).

Thus, the overall survival function is the mixture:

S(x) = p. f ₀ (x) + (1 -p). fi(x) (A.1 )

(f-i(x) may be different from survival function of the disease without a treatment because of side effects of the treatment).

We assumed that this relation is true for the medians (It's indeed true for the means but means are less usual in survival analysis, and simulations show that the bias for the mixture medians is not so substantial): OS = p.RS + (1 -p).NRS (A.2)

OS (Overall Survival) is the median time from the treatment to death/loss of follow up.

RS (Response Survival) is the median time from response to death/loss of follow up in the responders group.

NRS (Non Response Survival) is the median time from the treatment to death/loss of follow up in the non-responders group.

In studies, OS and CBR are frequently reported but not RS and NRS.

However there is another parameter which is frequently reported in

oncological studies and which may help reproduce RS and NRS (besides it being an important efficiency value by itself):

The PFS (Progression Free Survival) or TTP (Time To Progression) is the median time from the treatment to the progression of the disease.

In the best case the DCB (Duration of Clinical Benefit) is also reported, which is the PFS of the participants defined like CBR.

In the second best case the Duration of Response is reported, which is the PFS of the participants defined like ORR.

We will show below how we can graphically (and approximately) estimate the DCB from a Kaplan-Meier (KM) curve on PFS or Duration of Response.

We hereafter called the DCB: "RTD"

Let's assume:

RS = RTD + NRS (A.3)

(After the progression a responder is like a non-responder after the treatment)

Hence: OS = p.RS + (1 -p).NRS = p.(RS-NRS) + NRS

OS = p. RTD +NRS (A.4) OS - (l - p)x NRS _ OS - (l - p)x (OS - _P x RTD

p

RTD+ OS-p x RTD = OS + (l -p)x RTD (A.5)

We use this model twice. In a first time we estimate the general RS and NRS of the sample study with p = CBR. In a second time, if we get a probability to respond less naive than p, called here Molecular Prediction of Response or MPR, we can specifically estimate a more accurate OS for our patient (called here effOS):

i.e. After getting OS, RTD and CBR from the study, we successively compute:

RS = OS + (1- CBR).RTD NRS = RS - RTD effOS = MPR.RS + (l-MPR).NRS (A.6)

Let's describe how we estimate our RTD from Kaplan-Meier curve:

Without taking into account participants who die before progression, the RTD is the median in a KM curve that would ignore the first PD steps: For the graph of Fig. 7, the paper reports PFS =15m and CBR=73% . So if we remove the top 27% that are PD, the RTD is the median of the remained curve and we consider the survival time of CBR/2, here 0.73/2 = 0.365 gives 29m.

Now, the MPR (Molecular Prediction of Response) is based on molecular results of the tumor lab analysis like the following: At first, we attribute to the patient a "qualitative MPR": MPR.QI in {-1 , 0, 1 }

If we don't know about any molecular markers that can predict a clinical response, than: MPR.QI = 0

If the molecular analysis predicts no clinical response: MPR.QI= -1 , and if the molecular analysis predicts some clinical response (including stable disease) MPR.QI=1 .

MPR.QI = 0 => we use MPR=CBR

MPR.QI = -1 => we admit in a deterministic way that MPR = 0 + sd(CBR)

MPR.QI = 1 => we distinguish 3 cases. Before, let's define the events "the marker exists" Λ+ ,"the marker doesn't exist" Λ- and the events "the treatment works (including stable disease) " Ω+ ,"the treatment doesn't work" Ω-, and we assume (it's biologically defendable)

Λ- => Ω- (or Ω+ => Λ+) · In the best case if the article gives Ρ(Ω+|Λ+) we use it.

• If we have relevant knowledge on Ρ(Λ+) we use a Bayesian estimation

Ρ(Ω+|Λ+) = CBR/ Ρ(Λ+) (2.6)

• If we do not have either, we'll use a (heuristic) deterministic estimation

Ρ(Ω+|Λ+) = CBR + a.sd(RR) (2.7) where a=3 for chemotherapy treatments and 2 for biological treatments.

Summary for estimating the MPR.QI

parameter MPR :

1 0 -1

The article

Knowledge Ρ(Ω+|Λ+) CBR , Ρ(Ω+|Λ-) from articles gives Sd calculable sd(CBR) Sd calculable Ρ(Ω+|Λ+) other CBR/ Ρ(Λ+)

articles give

Sd calculable

Ρ(Λ+)

sd(CBR) no CBR+a.sd(CBR)

knowledge

ae{2,3} Without sd

Table 1

To estimate the standard deviations we use confidence level given in paper for median time and for proportions like CBR we use

The three AEV are finally: the RTD adjusted or ART=RTD.MPR, the effOS and a subjective OS (we choose explicitly in which group the patient is) i.e.

ARS=MPR.RS if MPR.QI=1 APS=MPR.NRS if MPR.QI=-1 .

(If MPR.QI=0, the third AEV is not used) Appendix II

In this appendix, we illustrate the model and the algorithm by a true example.

The treatments alternatives for a patient X and their efficacy estimators

reported in clinical trial studies on populations in which the patient could be eligible (according to the exclusion/inclusion criteria) are summarized in the following table:

Table 2 ^* from TTP or PFS and graphical KM considerations described in Appendix I

The standard errors of the estimators are restored from confidence level and sample sizes (For the CBR proportions the sample size is enough to compute its standard error)

The computation of the RS and NRS is easy from the CBR, TTP and OS

estimators by the assumptions RS = OS + (1 - CBR). RTD; NRS = RS - RTD explained in Appendix I and we get (in months):

Treatment

alternative RS NRS

Raltitrexed 11.48 6.68

Dacarbasine 9.12 7.21

S-1 11.65 8.85

Regorafenib 7.52 5.62

Cetuximab 13.18 6.58

Afatinib 16.07 12.54 Everolimus 7.12 5.32

Tivantinib 8.69 5.99

Table 3

The standard errors are given by the theoretical equations:

Var(NRS)= Var(OS - RTD)= Var(OS)+ Var(RTD) > a(NRS)= ^a ² (OS)+a ² (RTD)

Var[( - CBR)x RTD]= Var[CBR x RTD]= Var(ART)

Var(RS)= Var(OS + (l - CBR)x RTD)= Var(OS)+ Far((l - CBR)x RTD)= VarpS)+ Var(ART)= Var(NRS) The molecular predictions of the subject with the appropriated efficacy

estimators or prevalence estimator from clinico-molecular and epidemiologic studies give the whole information that we need to obtain the final MPR. See table 1 .

Table 4 From RTD (table 2), MPR (table 4), RS and NRS (table 3) we get our three AEV by ART = MPR. RTD; effOS = MPR.RS + (l-MPR).NRS; subjOS = MPR.RS if MPR.QI=1 or MPR.NRS if MPR.QI—1. Our three AEV are in the following table and constitutes the initialization of the 3 algorithms:

The subjOS is not available for Regorafenib because there is no molecular prediction for it.

Let's describe the first steps for one of the three algorithms which have to be done by the ART example:

Step 0

The alternatives (designated here by their references) are ordered by ascending ART, and CLR are computed

. the first one between treatments 5 and 4 is

ref n ART se(ART) CLR Average CLR

5 39 0.594 0.527 62.31

4 505 0.779 0.265 52.528

2 68 0.836 0.859 56.199

7 71 0.971 0.104 65.1 18 58.53986 8 71 1.35 0.97 65.289

3 27 2.182 1.881 54.088

1 20 2.4 0.986 54.247

6 41 2.516 0.459 NA

Step 1

The smallest CLR is between treatments described at ref #4 and #2 (see above) so we pooled them in a new ART with a new standard error and a new sample size.

A m- - 505.0.779 + 68.0.836 . .. . . .

The new ART is 0.785764 = and its new standard error is

505 + 68

Step 2

In the same way we're pooling here treatments 3 and 1 ref n ART se(ART) CLR Average CLR

5 39 0.594 0.527 61.535

4 and 2 573 0.785764 0.387037 67.803

7 71 0.971 0.104 65.118

8 71 1.35 0.97 69.203 63.9058 3 and 1 47 2.274766 1.567762 55.87

6 41 2.516 0.459 NA

Step 3

In the same way we're pooling here treatment group 3 and 1 with treatment 6

The last step (6) gives:

The step functions remaining from the entire procedure are depicted in Fig. 8.

The descending function represents the decreasing number of treatment groups and the ascending function represents the increasing average CLR. The x-axis is constituted by the CLR intervals of the steps and the unit of the y-axis is the interval (2-8) corresponding to the number of treatments and in which the average CLR was projected.

The maximal product is given by step 3 (the intersection in the graph above) The ART are projected on the grade interval (1 , 10) and we get:

Raltitrexed 10

Dacarbasine 2

S-1 10

Regorafenib 2 Cetuximab 1

Afatinib 10

Everolimus 2.9

Tivantinib 4.8

After the 3 algorithms are done we can summarize 3 grades for each alternative and the final grade is the weighted mean of them (The weights are 0.5 for ART, 0.45 for effOS and 0.05 for subjOS when available)

Treatment grade grade grade

alternative ART eff.OS subjOS Total. Efficiency. Score

Raltitrexed 10 3.8 3.2 6.7

Dacarbasine 2 2.2 1 1.6

S-1 10 5.8 7.1 7.8

Regorafenib 2 1 NA 1

Cetuximab 1 2.2 3.2 1.2

Afatinib 10 10 10 10

Everolimus 2.9 1 1 1.5

Tivantinib 4.8 2.2 1 3.1