Variable transformation based on exponential decay

[Alireza_Boloori]

Hello everyone,

I have a dataset of patients who're prescribed two medications at different time intervals. This is an example:

Patient ID    Drug1_dose    Drug2_dose    Time_start    Time_stop

1                  10                    5                     1/1/2017        1/5/2017

1                  0                      7                     2/17/2017     2/25/2017

2                  3                      0                     1/15/2017     2/10/2017

....

Now, I want to apply a weight based on an "exponential decay" function to the drug doses, so that they can be extended over an extra period. Here is the function: 

$exp{-(x^2)/(2a^2)}$

 where "x" is the number of days elapsed since when the drug was prescribed, and a is a constant (a=0.849322). For example, for patient 1, drug 1, the extended drug dose is calculated as follows (extended time windows = 3 days):

Date	1/1/17	1/2/17	1/3/17	1/4/17	1/5/17	1/6/17	1/7/17	1/8/17
Dose_original	10	10	10	10	10
Weight (from 1/1/17)	0	0.5	0.0625	0.001953
Weight (from 1/2/17)		0	0.5	0.0625	0.001953
Weight (from 1/3/17)			0	0.5	0.0625	0.001953
Weight (from 1/4/17)				0	0.5	0.0625	0.001953
Weight (from 1/5/17)					0	0.5	0.0625	0.001953
Total	10	15	15.625	15.64453	15.64453	5.644534	0.644532	0.019531

0.5 = exp{-(1^2)/(2*0.849322^2)}

0.0625 = exp{-(2^2)/(2*0.849322^2)}

0.001953= exp{-(3^2)/(2*0.849322^2)}

I was wondering if there is any way to make this transformation. I Know that the number of rows (observations) may get bigger hyper-exponentially, but I need to have this transformation. 

Thank you for any feedback/thoughts, in advance!

[mkeintz]

1. Just in your small sample, you have 56 implied dates (1/1 through 2/25).  Do you really want 56 columns (as is implied in your sample output)? More generally, if your total date range is 1 year, do you want 365 columns?  Perhaps you might prefer a data set organized with just 5 vars, and as many rows as needed to accomodate your decaying weight:
   
   

   data mydata;
    input id dose (from_date to_date) (:mmddyy10.) weight;
    format from_date to_date date9.;
   datalines;
   1 10 1/1/2017 1/1/2017 1
   1 10 1/1/2017 1/2/2017 .5
   1 10 1/1/2017 1/3/2017 .0625
   1 10 1/1/2017 1/4/2017 .001953
   1 10 1/2/2017 1/2/2017 1
   1 10 1/2/2017 1/3/2017 .5
   1 10 1/2/2017 1/4/2017 .0625
   1 10 1/2/2017 1/5/2017 .001953
   1 10 1/3/2017 1/3/2017 1
   1 10 1/3/2017 1/4/2017 .5
   1 10 1/3/2017 1/5/2017 .0625
   1 10 1/3/2017 1/6/2017 .001953
   1 10 1/4/2017 1/4/2017 1
   1 10 1/4/2017 1/5/2017 .5
   1 10 1/4/2017 1/6/2017 .0625
   1 10 1/4/2017 1/7/2017 .001953
   1 10 1/5/2017 1/5/2017 1
   1 10 1/5/2017 1/6/2017 .5
   1 10 1/5/2017 1/7/2017 .0625
   1 10 1/5/2017 1/8/2017 .001953
   run;
   

   You could then run a proc tabulate with DOSE as analysis var, class vars of from_date and to_date, and weight as the weighting vars.  It would yield a table laid out like your requested result.
   

   proc tabulate data=mydata;
    by id;
    class from_date to_date;
    var dose;
    weight weight;
    tables to_date all, from_date*sum*dose=' '*f=best10.6 /rts=12;
   run;

   
   If so, then making MYDATA from your original dataset is fairly simple.       
2. Your sample decay is censored after 3 days.  Is that your intention?

[Alireza_Boloori]

@mkeintz Thank you for the response! However, there are a few points that need to be mentioned:

1. The whole idea behind this step is to prepare the data for conducting "survival analysis."
2. So, although all 3 years are rquired, each patient was prescribed for a specific part of this 3-year period (some had multiple time intervals).
3. The second table I provided was for clarfiying my point (it was not the output I'm looking for). The output I'm looking for has the same structrue as my first table (with the difference that now I've extended dose over extended time period).
4. The "mydata" dataset you created at the beginning of your response cannot be helpful, because (i) it's not generated by a macro, and (ii) for a very large-size data (like mine with millions of patients), it'd be impossible to generate it.
5. The extended time windows for this example was 3 days. However, I may have to change this threshold to try other values. That's another reason that I need a macro.

[PGStats]

I'm puzzled. Your decay equation is not not one of exponential decay because it is quadratic in time. It is as if molecules that have been in the system for a longer period are eliminated faster than freshly administered drugs. How did you get that decay equation?

[Alireza_Boloori]

At any day a drug is prescribed, it'll have delayed exposure in later days. That's where the exponential decay function is applied (here is the source: https://www.ncbi.nlm.nih.gov/pubmed/16549262). The numbers you see (e.g., 15.625) is the accumulation of all delayed exposure attributed to the previous days (whose numbers are censored by the time window).

[PGStats]

Sorry, I don't have 35.95USD to spend on getting access to that article. Do the authors discuss the more than exponential decay rate equation?

I'm picking on this because truly exponential decay is much easier to handle. At any point in time, you only need to know the concentration at that time to predict the concentration at later times. With your equation, you need to know the full drug administration history to predict the concentration at later times. 

[Alireza_Boloori]

My bad! I thought the article would be available under an academic access. Regarding the full drug admin history, I provided that for the patient # 1 in my first table (please see above). I'm not able to provide the history for all patients, because as I said before, there are millions of patients in this data.

[PGStats · Accepted solution]

This is how to do this with an array:

data have;
input ID D1 D2 (start stop) (:mmddyy10.);
format start stop yymmdd10.;
datalines;
1 10 5 1/1/2017 1/5/2017
1 0 7 2/17/2017 2/25/2017
2 3 0 1/15/2017 2/10/2017
;

/* "decay" constant */
%let a = 0.84932180028801904272150283410289;
/* Decay window (days) */
%let dw = 10;
/* Define total date range, allowing for a decay period at the end */
%let ldate=01jan2017;
%let hdate=10jan2018;

%let lbound=%sysfunc(int("&ldate"d));
%let hbound=%sysfunc(int("&hdate"d));

data want;
format date yymmdd10.;
array dose(&lbound.:&hbound.) _temporary_;
call missing(of dose{*});
do until(last.id);
    set have; by id;
    do date = start to stop;
        do _d = 0 to &dw;
            dose{date + _d} + (D1 * exp(-_d**2/(2*&a**2)));
            end;
        end;
    end;
do date = lbound(dose) to hbound(dose);
    if dose{date} > 0.00001 then do;
        drug_dose = dose{date};
        output;
        end;
    end;
keep id date drug_dose;
run;

proc print; by id; id id; run;