Electrónica Nova Scientia Radioterapia y el papel del sistema inmunitario Radiotherapy and the role of the immune system

Using a population dynamics formulation to model radiotherapy treatments a large amount of patients under treatment was mimicked. From the obtained results a method to change the number of radiation sessions while keeping the optimized treatment is found. Also a power law relating immune system tumor effciency ratio (ISTER) for a given patient and maximum success probability of radiotherapy treatment is obtained.


Introduction
The relevance of immune system tumor interaction in a radiotherapy treatment (RT) and the importance of its study have been highlighted in a previous work [1].The role of immune system in tumor control has been described in [2] and widely discussed in further works [3][4][5][6][7][8].However, in our opinion, the need of its study and characterization to plan a treatment has not been fully understood.
Mathematical models are widely used to describe tumoral behaviour withor without treatment.
Such models have been used to explain and to describe immune depression effects [9], the role of time delays in immune response [10,11], radiovirotherapy treatments [12], chemotherapy treatments [13] and even surgery [14].
There is an increasing concern about finding the suitable planning that maximizes the outcome of a radiotherapy (RT) treatment [15].Generally speaking, the problem of how to apply radiation treatments in an optimal way taking into account the influence of external factors is, in our opinion, far from being fully solved, specially in what concerns the derivation of general decision tools to organize and adopt radiotherapy protocols in actual clinical conditions.
In the present work we investigate how to apply the optimum dosage and to adapt this optimal treatment to other protocols with a different number of sessions.This could be useful to oncology services to plan the RT treatments without decreasing treatment performance.Also the oncologists could adjust the radiation dosage in order to minimize the damage caused to the surrounding tissue.

Model
We will use a Lotka-Volterra like model to describe the tumor evolution based on some assumptions.Tumor cells growth X (as usual, a dot over a quantity represents its time derivative) depends on the current tumor population as aX and its interaction with lymphocytes, -bXY .Lymphocyte population grows due to tumor-immune system interaction, dXY , and contribution to exponential decay, -fY , due to natural cell death.The tumor is assumed to secrete interleukin which produces an immunity depression effect [16,17], -kX.In this model a constant flow, u, of lymphocytes arrives from the immune system.So, we model tumor-immune system interaction using the already known equations [9]: The inclusion of both populations affection brings a new equation for the tumor non clonogenic cells [12], Z, originated from radiation damage.We assume that those cells will also stimulate the lymphocytes population, as pZY , will decay exponentially as -rZ due to the death of damaged cells, and as -qZY due to the interaction with immune system.Finally we arrive to where ( ) ( )


. All parameters can be estimated and interpreted by a similar procedure as in [11].
A linear stability analysis of the system (3) shows that tumor will vanish to 0 L = (0; ; 0) Figure 1 shows stable and unstable regions of Eqs. ( 3) and highlights region III on which this work will focus.If the system is L 0 -stable and initial tumour size is small enough, then the radiation treatment is unnecessary, whereas if tumour size is large enough, then the treatment will take it closer to L 0 .
The L 1 controlled growth state will be reached only if both parameters fulfill the same condition, in other words, if   and  are both greater or smaller than unity at the same time.Any other condition makes the first component of L 1 < 0, and even when the stable point mathematically exists, it can not be approximated from realistic initial conditions (that should remain positive along the simulation time).For those patients with 1   and / < 1  , the main effects of the tumor will be the depression of immune system, they will perform badly according to Karnofsky performance scale [18] and will not fulfill physical requirements to be subject to treatment.

However for 1
   , tumor will grow and tumor eradication will be achieved only by bringing the system close enough to L 0 , so that the immune system can get rid of the tumor.
The chosen characteristic time and the dimensionless parameters allow us to give a very intuitive interpretation of the critical parameters of Eqs.(3).We can see /  as the effciency of immune system over tumor growth and  as the "deficiency" of the immune system due to tumor growth.
It is also easy to see that radiation treatments do not change the stability conditions of the system (2), since radiotherapy does not change tumor or lymphocytes growth rate, but can drive the number of both kind of cells to very small values.This means that, for the chosen region in the parameter space, any remaining tumor cells will eventually reproduce and grow after the end of the treatment.
Although Eqs. (3) allow for infinitesimal x values, in real systems when the number of tumor cells becomes small enough, immune system may kill them.In other cases, when a few tumor cells survive, they can cause tumor regrowth.It is known that this behavior is almost independent on tumor size [19].As an estimation we will assume that if the immune system has a high efficiency it will have higher probabilities of kill the tumor remaining cells and if the immune system has a low efficiency the probability of kill those tumor cells will be low.The simplest linear approximation to deal with this effect is to assume: as the probability of tumor regression.If no regression occurs, tumor will regrow.We must realize that when the immunodepression term -kX of system (1) is larger than the flux of lymphocytes u, an unattended method to solve the equations will fail to describe a biological system [3].However, this is a simple way to introduce the immunodepression effect in our equations and can be considered as the first order Taylor approximation for a more general nonlinear function.Furthermore, it provides us a very simple portrait for the parameter space that allows to select the target tumors of radiotherapy treatments.Finally, the chance of lymphocyte population to become zero gives us a natural cutoff for our system.Then, whenever lymphocyte population becomes zero we will assume the tumor escapes lymphocyte control and grows limited only by space and nutrient considerations.At this point, we consider treatment has failed.

Treatment optimization
A reference protocol was simulated in [1] while coefficients entering equation 3 were varied at random among admissible values [11].The simulation, then, covers a wide range of tumors and supply general results useful to clinical practice.
One million of "virtual patients" under treatment were simulated.Defining the probability of treatment success (P s ) as the fraction of "virtual patients" with no tumor at the end of treatment, it could be represented as a function of tumor cell survival factor and the Immune System -Tumor Effciency Rate [1], (ISTER), defined as ISTER = /  .
As shown in [1], the long term survival of patients will not improve with higher doses of radiation, on the contrary, it is possible to get the maximum success probability for our reference protocol at intermediate doses.Our model predicts that higher doses do not improve therapy outcomes.

Simulation
We can mimic different radiation treatments with Eqs.(3) to simulate tumor evolution.To follow radiotherapy treatment in a realistic way, we apply a radiation session every workday (5 consecutive days) and none in weekends (2 days).All treatments [20,15] begin the tenth day, and apply a variable number of sessions, N, from 20 to 40 for each patient.
We could take up a sort of tumors, i.e. breast or colon, and use an experimental expression to calculate the survival fraction of cells for a 2 Gy dosage.Instead, we preferred to proceed in a more general way, and ignoring the radiation dose, we have taken a random value for the survival fractions of each patient tumor.We have generated several virtual patients under treatment taking different values for the parameter values in Eqs.(3) and use a fourth order Runge-Kutta method [21] to integrate them.
To reproduce tumor evolution resembling that of a clinical case, one needs to calculate the correct values of the cofficients appearing in Eqs.(3).Numerical estimation of these coefficients was already made in [11] (and also in [8] for a slightly different model), based on clinically available data, showing a possible procedure for clinical professionals to estimate their values.
A statistical study of the dependence of treatment success on the dosage and number of radiation sessions was performed.Due to the wide range of possible parameter values in Eqs.(3), their values have been drawn randomly from a log-normal distribution, to avoid negative values, but keeping the immune system efficiency ( /  ) always smaller than 1.Survival factors [19, 22]   are also taken as random values.As initial conditions we assumed, for simplicity, that the number of tumor cells is higher than the number of lymphocytes and that both initial populations are distributed as normal random numbers.

Results and discussion
To describe the success probability for a given treatment let us define the effect potential for tumor cells as ln( ) where t S is, as before, the survival fraction of cells.This quantity grows monotonously with the radiation dose.If an exponential model like the linear or the LQ model is used to describe the interaction of radiation with living tissues, then the potential becomes the tissue effect.
Otherwise, like in the single-hit multi-target inactivation model, it has not a direct interpretation but remains useful to describe the radiation damage.
The treatment success probability, P s , defined as the fraction of patients with no tumor at the end of treatment, is represented as a function of ISTER, N and  .
As expected, for each value of N we obtained similar results as obtained in [1].Furthermore, when a fixed value of the ISTER is taken, the surface represented in figure 2 is obtained.This shows that whenever the value of N increases, the optimized value of P S can be obtained with a lower value of  per session.This allows to find the exact values of  that optimize P S , for each value of N, as represented in figure 3.
where R is constant for a fixed ISTER value and max P is the maximum success probability the treatment is able to get, given  and N values.Those R values then could be fitted as a power law of the ISTER, Our simulations allow to obtain m = 1,1 0,1  and A = 17  1 bringing a method to obtain an approximate value of the maximum success probability that a radiotherapy treatment could reach if the approximate ISTER value is known.The parameter R allows to compare without further information the optimal treatments that could be applied.Tumors in zone 1 have a high value of R and even with high radiation doses will not get a good success probability.Taking a standard protocol with 30 sessions as a reference, figure 5 shows that each radiation session must get a value of  = 0,55 and then according to expression 7 a treatment will get less than 14% of success probability.However an R value between 20 and 30 like in zone 2, guarantee that the final success probability will be high if an optimized treatment is applied.The same standard protocol of 30 sessions must apply a radiation dosage equivalent to  = 0,6 to get the maximum success probability.Here R = 24 guarantee that if an optimum treatment is applied, a success probability of 75% could be achieved.
Thus, large values of R mean that the success probability of a treatment is low in comparison with the same treatment corresponding to small values of R. The main factor in this analysis is, then, the ISTER.This shows how strongly the immune system influences the results of the radiation treatments and illustrate the importance of considering it as a main factor for radiotherapy.

X
= aX -bXY Y = dXY -fY -kX + u (1) as a result of RT treatment a fraction of cells of both populations, lymphocytes and tumor cells are affected.The lymphocytes population get quickly reduced in a fraction l B of cells in what is called the interphase death.However, affected tumor cells will not die immediately but a fraction t B loses its reproductive capacity and dies in the mitosis process.The fraction of affected cells in both populations is the complement of the survival factors t of tumor cells and lymphocytes affected by radiation per unit time.T n are the time instants when radiation doses are applied and (T -Tn)  denotes Dirac's delta function centered at Tn.A dimensionless system can be easily obtained taking the tumor duplication time 1 influences) as the characteristic time, so we introduce the dimensionless time substitutions X = ax/d, Y = ay/b, Z = az/d, we get t t x = x -xy -( )x y = xy + zy -y -(x + z) + -( )y z

Figure 2 :
Figure 2: P S representation for ISTER = 0,7 as a function of the effect potential and the number of radiation sessions.Closer to black means lower, yellow means closer to 0,7

Figure 3 :
Figure 3: Effect potential per session that maximize the success probability Ps for each value of N plotted for some values of ISTER parameter.

Figure 4 :
Figure 4: Schematic representation of zones grouping example tumors.

Figure 5 :
Figure 5: Effect potential per session that maximize the success probability P S for each value of N plotted for the zones 1 and 2 of parameter space showed in figure 4.