7 replies to this topic

[Lance Wollenberg]

Is there a way to elegantly code an expression to do a time above threshold where a concentration is above a set value in PML?

 

I know how to do this with brute force (Simulate, set an if-else statement in R, summate by patient), but this is not ideal.

 

Thanks!

Lance 

[serge guzy]

Dear Lance

I am attaching a project that answers your question.

 

The code is for a simple 1 compartment model as below:

The idea is that you must integrate the time from 0 to the final time point only when the predicted concentration is above a threshold value.

 

In order to do that you first need to define an index that will take the value of 1 if the predicted concentration above the threshold and 0 otherwise.

Then you write

deriv(tthreshold=index)

 

tthresholdis the time above the threshold level and if you were always above the threshold index would be 1 all the time and derive(tthreshold=1) would be equal to the actual time at each time.

 

This comes from the definition of a derivative.

 

dX/dt=1   dx=dt which means that the change of your response (x) equal the  change in time.

however when  C is below the threshold, index=0 and therefore during that time the threshold does not accumulate.

the overall result is that for each patient you can follow at each time how much time they were above the threshold.

Obviously this time can never be larger than the current time. BEST SCENAIRO IS EQUAL WHEN YOUA RE ALWASY ABOVE THE THRESHOLD.

I am attaching an example project to help you.

best Regards

Serge

 

 

the code for a compt model  IV input

 

test(){
 deriv(A1 = - (Cl * C))
 urinecpt(A0 = (Cl * C))
 C = A1 / V
threshold=60

index=(C>threshold?1:0)
 deriv(tthreshold=index)
 tthresholdgraph=tthreshold
 dosepoint(A1, idosevar = A1Dose, infdosevar = A1InfDose, infratevar = A1InfRate)
 error(CEps = 0.1)
 observe(CObs = C * (1 + CEps))
 stparm(V = tvV * exp(nV))
 stparm(Cl = tvCl * exp(nCl))
 fixef(tvV = c(, 1, ))
 fixef(tvCl = c(, 0.2, ))
 ranef(diag(nV, nCl) = c(0.1, 0.1))
}

Attached Files

•  time_above_threshold.phxproj   911.03KB   372 downloads

[Lance Wollenberg]

Thanks for the always helpful suggestions Serge!

 

Lance

[serge guzy]

My pleasure as always.

Serge

[Simon Davis]

Lance, I read this request  as that you wanted to know time spent above  e.g Cp of 60

 

In which case you could code it as one line.

 

deriv(t_above60 = (C>60)*t)

 

Anf then add it to your existing table statement that would give you cumulative time spent above the threshold conc of e.g. 60 in Table01

 

or You could ask for this data at e.g. 999 h in own table e.g Table 02.

 

I'm not sure if this is clearer for you/others than Serge's code - or maybe I have got the wrong end of the stick ?

 

 Simon  time_above_threshold_saad.phxproj   920.71KB   373 downloads

 

Attached Thumbnails

• 

[Lance Wollenberg]

Thanks Simon!

[Lance Wollenberg]

Hi Simon,

 

I actually think the correct equation should be deriv(t_above60 = (C>60)), so this would give you the amount of time where C>60. 

 

deriv(t_above60 = (C>60) * t)

 

In the equation that you proposed deriv(t_above60 = (C>60) * t), this would give you the AUC of the Time above the threshold. Agreed?

 

Thanks!

Lance

[Simon Davis]

Sorry Lance,you were absolutely correct about the time_above_t however the AUC above this threshold only requires a few more steps. here is some annotated code to cover most eventualities - should all be correct now  ;0)

 

threshold=60            # set a threshold conc
index=(C>threshold?1:0) # set index to 1 when C greater than threshold (i.e. 60)
    
deriv(AUCinf = C)             # simple AUCinf                     
deriv(AUCt6 = (t>0)*(t<6)*C)  # calc AUC from time 0 to time 6    
deriv(t_above60 = (C>60))     # calc time greater than 60         

# step 1: define AUC above 60 in the whole time range
# this is AUC cumulated over time (from 0 to current time) when C >60
# note it is ALL the AUC
deriv(auc_above60all = (C>60)*C)

# now we use previously defined the overall time above 60
# and multiply t_above60 by the threshold value i.e. 60
    rectangle=threshold*t_above60
# now we substract AUCABOVE60all  by this rectangle to give us the AUCabove60 only
auc_above60=auc_above60all-rectangle          

 

that would be the pink area below, Simon

 

Attached Thumbnails

•