Calculate I(d) after the slope of an IDER becomes 0

[eghuang]

From Ray's notes, 07-17-2017:

"The figure shows computations for a 50-50 mixture of a high LET Fe beam and a low LET 1H1 or 2He4 beam. The IDER for the Fe beam is the IDER we have developed and actually plan to use because it is a big improvement over the previous models. It is calibrated with our actual data. It never gets bigger than 1 because the endpoint “fraction of mice that have at least one HG tumor” cannot be greater than 1. That is different from our previous figure but not in itself a problem because of course we are not interested in mixture effects bigger than 1 either.

The low LET light ion IDER is a toy IDER calibrated with our actual data and using it requires care. For some reason the light GCR ion data strongly indicate an upper limit of about 0.6-0.8 for this endpoint. Our toy IDER for this data has a maximum of 0.8 which it never exceeds. So once the HZE ion drives I(d) above 0.8 uniroot can not find the now non-existent intersection point and I(d) becomes undefined, even though we are interested in effects between .8 and 1. That kind of problem occurs often and causes trouble.

One solution is to argue as follows. As we approach I=0.8 from below ion 2 is already doing very little. At 0.78 uniroot can find the intersection point and we can calculate the slope m2 there to add 0.5m2 to dI/dd.But m2 will obviously be very small since the curve is so flat there (whereas at I=0.8 slope m1 is still very large). So we simply extend our definition to say that at I=0.8 or more ion 2 does nothing at all; all further changes to I(d) come from ion 1, which eventually drives I(d) up to 1, but slowly because ion2 is only acting at half strength.

So one of the two main steps remaining for our whole project is to program this idea carefully using actual instead of toy IDER for the low LET contributions. Care is required. If we let ion 2 stop at I=0.78 that is a bit too soon. If we try to go to I=0.79999 uniroot may miss by mistake and think there is no solution. But that is just a technical programming issue. In principle we have a solution and only need to implement it.

This solution we will use is not the best solution. If ion 2 hates effects above 0.8 why doesn’t it put up a fight when ion 1 drives I(d) above 0.8 and force a compromise upper limit between 0.8 and 1? There is a way to do that but we will not use it for our first paper."

[rainersachs]

Edward and Mark: it is time to address this issue. Have you guys been working on it? Do you need a meeting on it?

[eghuang]

I wanted to clarify a few things:

1. What kind of script do we want? Are we looking for a general function that takes in two IDERS, a dataframe, and outputs a MIXDER based on your idea? Are we just writing several lines of code specific to our model and current dataset? What data points are we using to make these calculations?
2. If we are writing a function, how do we want to specify the stopping point for ion x?

Generally, additional information on the inputs and outputs would be extremely helpful in implementing this idea.

[rainersachs]

The answers to question 1 are the following: (a) Any script that does the job, whether based on R functions, loops, or anything else is OK. You should do whatever you think is best for R and see that it works. (b) We want it for mixtures containing as many as 8 components, but please start with 2-component mixtures. (c) Start with code specific to our model; other models can come later for later papers. (d) We are using our data; specifically the data frame, already in place that contains all the non-zero dose results and omits one of Fe data sets.

The answer to 2 is that there are two possible methods. (a) calculate the derivative analytically. the stopping point is just very slightly before the derivative becomes zero. (b calculate the derivative numerically. The stopping point is again the same.

Whichever works better but my guess is that will be (a).

IIRC I estimated the stopping point crudely with both (a) and (b) and convinced myself they were roughly the same answer. See some of the commented out parts of the script

[rainersachs]

I didn't mean to close the issue

[eghuang]

Okay, so whatever we write should identify the greater IDER (assuming we're only looking at two IDERs for now), the dosage at which we decide to exclude its effect on I(d), and its slope at that point. Does that sound right?

Also, How would I get the function for each IDER?

[rainersachs]

it should identify the lesser of the two IDERs, i.e. the one that never makes it up beyond about 60%, i.e the low LET IDER for our data set. The dosage at which we decide to exclude exclude its effect on the I(d) is always the dose at which the slope is (very slightly more) than zero. So the slope at that dose is (very slightly more) than zero. The function for each IDER is the function called our model, or our IDER, in the script. For the high LET ions it depends on dose and on the auxiliary explanatory variable L=LET. For the low LET ions [in general,fast ions for which Z<4; in our data set Z=1 (Hydrogen) or Z=2 (2He4, i.e. alpha particles)]. The script DIDER2 (2).R DIDER2 does this calculation for a toy example with simulated data and produces the graph you gave at the very top of this issue

[rainersachs]

I just added a script where Siranart solved our problem but with now-obsolete IDER. Also Dae solved our problem for a toy case in DIDER2(2).R If you have a solution that works for N=2 you probably should ignore their methods; but otherwise you may find their ideas useful.

[eghuang]

A script has been written that satisfies the specifications of this issue.