sachsURAP / LSSR-2019

Data and code for the Life Sciences in Space Research 2019 paper. Concerns modeling murine Harderian gland tumorigenesis induced by mixed radiation fields.
GNU General Public License v3.0
1 stars 1 forks source link

Other Projects and Topics of Interest #6

Open eghuang opened 6 years ago

eghuang commented 6 years ago

From Ray:

This is a new issue for Mark to work on as his project for the rest of the semester after he finishes off the Lam report. The project is to study the generalized IEA equation using an approach based on the theory of functions of one complex variable. Yimin should make it one of his two main projects, the other being on call to help Ray and Edward with quality control on the master HG tumorigenesis script. Claire and Edward should participate iff time permits.

Claire has uncovered a new approach that I think might (or might not) lead to an actual break through. Noticing similarities between defining one-agent DERs via an autonomous ODE initial value problem and the R package that can integrate functions (and differentiate them) she suggests that instead of using dE/dd=F(E), and using dI/dd=Sum== Sum(from k=1 to the number of components of a mixture) F_k(E), one can simply use integrate on dd/dE=1/F(e), and on dd/dI=1/Sum. I'm pretty sure that at least in a sufficiently small neighborhood of the origin Claire's approach and the one we had started using, equivalent even for the oddball case F(E=0)=0.

However if you ask a modeler to guess, on the basis of experimental information plus biophysical knowledge or intuition, a suitable equation for 1/F(E) rather than guess a suitable equation for F(E) the modeler may well make different guesses. For example suppose the modeler starts with a finite polynomial near the origin and then decides to generalize to a convergent polynomial with an infinite number of terms. The result is trivial for our old method: merely the Taylor series for our function which we took to be always holomorphic near the origin. But I suspect you could find examples where taking a similar limit on [1/(finite polynomial)] gives something different from 1/(Taylor series). More generally there may be cases where 1/[limit F(e)] != limit[1/F(E)]. Or there may be other situations. Unless I am overlooking some theorem to the contrary, Claire's approach is not, for a modeler, equivalent to our current approach. It may then be better, or worse, or just a new alternative. We should try to study this using lots of examples and counter-examples and/or asking colleagues.

Edward: can you arrange for Yimin, Mark, Claire, and me to get write permission and to get notifications as regards this new Github-issue? Mark and Yimin please comment or ask questions if you have time.

Thanks! Ray

rainersachs commented 6 years ago

Yimin: please read, correct, extend, comment on, push to this issues thread, the following comments. Mark and Claire ditto if you have time.

In our meeting yesterday, Yimin had some good ideas and very interesting results on the functions of a complex variable project. He convinced me that we may be able to confine attention to polynomials P that have only a finite number of terms:

(1) P= a_0+a_1z+a_2z^2+...+ a_N*z^N. N < infinity.

That would certainly be an enormous simplification and in my opinion would show that the entire project can most probably be carried through to a useful conclusion.

First, Yimin argued that since we are never interested in very large doses we should always confine our attention to a dose interval [0,A], where A >0 is finite. Correspondingly when we consider complex effects we should consider a compact subset of the complex-effect plane for z=E+iy, where E and y are real and i=sqrt(-1). On the basis of some examples I suspect the appropriate subset is a disk of radius A centered at the origin [i.e. complex effect==z = 0 +0*sqrt(-1)].

Second Yimin pointed out that any function analytic on such a disk can (because of compactness) be approximated uniformly and absolutely by a finite polynomial P (Eq. 1) to any desired accuracy if we choose N large enough and determine the coefficients a_n appropriately. Since we are interested in data which is itself not infinitely accurate there will always be a sufficiently small inaccuracy for any given data set.

Next he wrote a script that picks finite polynomials "at random". I fear there must be many inequivalent ways to define "at random" here and some definitions are probably much better for our purposes than others. Even so this approach is certainly a useful supplement to, and is probably better than, my previous method of trying to choose "interesting" polynomials and investigate their properties.

Then he produced a number of fascinating examples which in my opinion suggest conjectures that if correct already solve about half the criteria needed to insure our final goal: having enough generality to ensure we can reasonably model any data set (of scalar effects dependent on dose) that is likely to arise in practice, but having enough restrictions on our polynomials to allow us to make relevant conjectures about and perhaps find actual relevant theorems that hold for every allowed polynomial. I will give some more details later, but have to break off now, and will push these comments now to prevent GitHub from eating them.

yiminllin commented 6 years ago

The main idea are outlined above. I will write some scripts for random polynomial generation, interpolation, and a custom ODE solver for you to play around with in next few weeks.

rainersachs commented 6 years ago

Hi Yimin:

In our discussion did you mean that you used random polynomial functions of dose, or random polynomial functions of effect E? Your clever argument that we can restrict attention to polynomials instead of having to consider more general holomorphic functions is applicable to either case. So is the excellent idea of using random polynomials. But I think we have to consider complex effect space z=E +iy (with E and y real and i=sqrt(-1)), not complex dose space z=d+iy. Please read the .pdf file here to see how and why this can probably be made to work. One may have to replace the basic polynomial assumption with Claire's suggestion that instead one should assume F(E) is the reciprocal of a polynomial (1/E is a polynomial) or use both the .pdf approach and Claire's approach simultaneously. instructions_to_Yimin_re_AC_one_agent_DERs_3.20.2018.pdf

yiminllin commented 6 years ago

I generate random polynomial F_j(I) (A5.2 in the pdf you attached), which is based on the effect. It is random in the sense that the coefficient are chosen randomly in real numbers (but could have complex roots). So the effect space is real number in this case. I could easily extend it to generate polynomial with complex coefficient though. The dose space is real, because we are solving the ODE from d=0 to =1.

rainersachs commented 6 years ago

Excellent! I look forward to seeing more examples.

Continuing the summary of our discussion last Thursday and looking forward to future calculations, here are suggested items for discussion and correction by mouse pod students.

Yimin has generated some random N=4 polynomials, of the Eq. (1) form,

P= a_0+a_1z+a_2z^2+...+ a_N*z^N. N < infinity, z=E+iy.

We looked at the resulting one-agent DERs, and at some IEA I(d) no-synergy/antagonism baseline mixture DERs. We discussed some results and conjectures. The 1-agent results and conjectures for the case a_0 >0 are the following. (2a) One can prove the 1-ion DERs are positive, with positive slope dE/dd on some half-closed E interval [0,A), where A>0 and A = infinity is allowed. Here and in the following we are considering the real (i.e. x) axis in the complex-dose plane effect==z=x+iy, where i=sqrt(-1) and x=E. (2b) Henceforth, A will be chosen to be the largest value for which (2a) is true. (2c) One can prove that there are three and only three qualitatively different cases: (2a1) A is finite and effect approaches a finite positive value as dose approaches infinity (reasonable behavior for a one-agent DER). (2a2) A is finite and effect approaches +infinity as dose approaches A (should not, as far as I can see, be allowed in any one-agent DER modeling actual data). (2a3) A is +infinity and effect approaches infinity as dose approaches infinity (sometimes a useful idealization for experimental data sets that describe one-agent DERs and for some or all mixture components in data sets that describe mixed-agent action). (2d) One can prove that if the zero of P that is closest to the origin E=0=y (i.e. z=0) occurs for positive E and y=0 (i.e. on the positive real axis), then case (2a1) must hold. (Probably this is the case that occurs most frequently in the literature on modeling experimental dose-effect data.) (2e) One can prove that if one of the zeros of P that is closest to the origin E=0=y of the complex-effect plane lies on the imaginary axis (the vertical, y, axis). Then the unacceptable behavior (2a2) occurs. (2f) Cases where a zero that is closest to the complex-effect plane origin is neither pure real nor pure imaginary have not yet been studied. Perhaps the key question is whether the real part is positive or negative. (2g) Importantly, with N=4 we saw cases where the number of inflection points of the 1-agent DER is any non-negative integer less than 4. This strongly suggests that by choosing N reasonably large one can model all the various shapes that inspection of actual data and biomathematical reasoning (or mathematical-physics/chemistry reasoning) suggests.

Turning to mixtures of a finite number M>1 of different agents, suppose each mixture component is described by a polynomial P as above, where the values of N and of a_k can be different for different components. Importantly, we allow mixtures where some a_0 coefficients are positive, some are zero, and some are negative. Define SIGMA as (SUM from j=1 to j=M of a_0). Assume for brevity SIGMA is non-zero. We also assume each component is allowed, e.g. is not case (2a2) when a_0 for that particular component is positive. Then we are pretty sure the following beautiful result holds. Incremental effect additivity is well defined and smooth for an interval [0, B) of the total mixture effect E where B is non-zero and real, with B = + infinity or - infinity allowed. Probably |B| can be chosen as large as the distance sqrt(E^2+y^2)==z(complex conjugate of z) == |z| from the origin to the point on the real effect axis where (SUM from 1 to M)(r_jP_j(E)) is zero Examples show all sorts of possible behaviors, importantly including many cases where some of the components' DERs have slope sign opposite to the a_0 of that component. The intuitive explanation then is that other components have pulled the effect beyond what the opposite-slope-sign component could ever achieve if acting on its own, and at such large-in-absolute-magnitude effects the component acts opposite to the way it acts at small effects (e.g. if it is an inhibitor near dose=0=effect, it becomes an effector when effect has a large absolute value, or vice versa).

This kind of mixture baseline behavior is far more general and potentially useful than for any previously published synergy theory dealing with scalar effects. It should probably be combined with or replaced by Claire’s modification where 1/E is a finite polynomial rather than E being a finite polynomial. No examples of Claire’s modification have yet been studied but that should eventually be done.

rainersachs commented 6 years ago

To discuss this Thursday with Yimin.

  1. In Eq. (1) above, all coefficients a_k have to be real. Let's assume temporarily that a_0 is non-zero. Then we can assume without essential loss of generality that a_0=1. I suggest when analyzing examples you just set a_0=1 and choose the other coefficients at random.
  2. What exactly does at random mean. For example are the coefficients assumed independent with each Gaussian-distributed with a zero mean and a variance chosen in accordance with some rule. If so, what is the rule?
  3. In case N in equation (1) is N=2 one could obviously compare the random polynomial results with the complete classification for a_0=1 -- 2 complex conjugate roots, 2 distinct real roots of opposite sign, two distinct positive roots, a double root = +-1.
  4. In case N=3 a similar classification exists. I'm told that follows from Galois theory but don't really know what Galois theory is. The classification is very messy and can hardly be useful but maybe the Galois theory is relevant to our questions.
  5. Anyhow an enormous amount is known about polynomials, including not only the fact that at least one root exists but all sorts of other stuff. For example if the coefficients alternate in sign and N is even the number of positive roots is even; if they alternate and N is odd the number of negative roots is odd.
  6. If we use Claire's approach dE/dd = 1/P(E), where P(E) is a polynomial with N=1 and with a_0=1, we can integrate in the form d=(integral from zero to E)(P(E')dE')=Q(E), where Q is a polynomial with real coefficients, a_0=0, and N'=N+1. Now if we complexify E to E+iy, dose d suddenly becomes complex instead of real. How did we end up with a complex dose and is there an interpretation?
  7. In 6, are we dealing with the differential equation dd/dE=P(E), where E is a complex variable? A lot is known about such equations.
rainersachs commented 6 years ago

Hi Yimin:

To be clear, the main priority should still be helping Edward with the mouse script if he needs help, not with this complex variable project. However, I did think a little more about the complex variable and think I found a simplification.

As you said, getting the zeros of a polynomial from from its coefficients involves instabilities and ill-posed problems (e.g. Wikipedia "Wilkinson's polynomial"). But our questions concern the roots, not the coefficients: what pattern of roots goes with specific DER shapes? So we should pick the roots at random, not the coefficients. That makes for much greater stability, addresses our question much more directly, and leads to simpler scripts.

This could be implemented as follows.

  1. Choose a degree N. If N is even 1.1. Choose E and y independently from uniform [-u,u]. 1.2. Take two roots to be E +- iy, giving the factor E^2 + y^2 in the polynomial. For later use store this factor E^2+y^2, which specifies how far the root is from the origin in the complex plane. To avoid platform dependent instabilities discard cases where E^2+y^2 <10^-13 (say). 1.3 Iterate step 1.2 N/2 times. 1.4 Choose the polynomial to be A*[(PI from j=1 to N/2) (E_j^2+y_j^2)], where A is that unique normalization constant that makes a_0 for the polynomial = +1
  2. When N is odd, carry out steps 1.1-1.3, add step 2.4 by choosing the Nth root to be a random real number E_j, and then calculate the normalization factor A.
  3. Now integrate to find the 1-agent DER. Investigate how the qualitative shapes depend on the root locations. Can we get pretty much any shape by choosing N and root locations carefully? Is it true that roots closer to the origin are more important than roots further away? Is it true that there are are only 3 main types of qualitative shapes? Is it true that we can make the inflection points follow any predetermined pattern? Etc.