Coding INTEGER Function in NONMEM

Hi, I would like to build a model to describe a circadian rhythm. Basically I assume that the release rate is linearly increasing from Tmin to Tmax, and linearly decreasing from Tmax to Tmin in 24 hr period. My code is as follows (I tried to put integer operation in $DES, to convert every time point after 24hr to 0-24hr period). My run failed at the covariance step, and nonmem gave me a warning: (WARNING 68) THE INTEGER FUNCTION IS BEING USED OUTSIDE OF A SIMULATION BLOCK. IF THE INTEGER VALUE AFFECTS THE VALUE OF THE OBJECTIVE FUNCTION, THEN AN ERROR WILL PROBABLY OCCUR. How the warning will affect my nonmem output? Why the run failed? Any suggestions on the coding will be highly appreciated. Thank you. Best wishes, lili NONMEM CODE: ;PD Model KOUT = 0.4 ;Elimination rate constant (1/hr) RMAX = THETA(1)*EXP(ETA(1)) ;Maximum release rate (ng/ml/hr) TMIN = THETA(2)*EXP(ETA(2)) ;Time at minimum release rate (hr) TGAP = THETA(3)*EXP(ETA(3)) ;Time from minimum rate to maximum release rate (hr) BSL = THETA(4)*EXP(ETA(4)) ;Concentration at time zero (ng/ml) SD = THETA(5) ;CV% of residual TMAX = TMIN+TGAP A_0(1) = BSL $DES T1= T-24*INT(T/24) ;Convert time to 0-24 period, the MOD function is not working in NONMEM IF (T1.LE.TMIN) THEN ;Make the time start from Tmin T2=T1+24 ELSE T2=T1 ENDIF IF (T2.LE.TMAX) THEN KIN=RMAX*(T2-TMIN)/TGAP ELSE KIN=RMAX*(TMIN+24-T2)/(24-TGAP) ENDIF DADT(1)= KIN-KOUT*A(1) $ERROR IPRED=F Y=F*(1+SD*EPS(1)) IRES=DV-IPRED IWRES=IRES/SD/F

Lili, you create a discontinuity in the first derivative and that can be difficult to handle numerically. Other than that, nonmem's failure to execute the covariance step is not necessarily a sign of anything. I assume that your two halves, Tmin and Tmax, are not 12 hours apart and that is why you do not use a standard trigonometric function (sine or cosine) for modeling the circadian rhythm. Here is an alternative: Tweak the standard cosine function to split into two parts of different lengths. Basically you stretch the x-axis differently in the two parts. L2=24-L1 ; L1=length of first period, L2=length of second period T24 = (TIME-CSHIFT) - INT((TIME-CSHIFT)/24); modulo 24 T1=T24/L1*12 ; part 1: map from (0, L1) to (0, 12) T2=(T24-L1)/L2*12+12 ; part 2: map from (L1, L2) to (12, 24) IF(T24.LT.L1)THEN TX=T1 ELSE TX=T2 ENDIF ; step-function but at L1, T1=T2=12 and gradient=0 PI=3.141592653 CIRC = CAMP*COS(TX/24*2*PI) Now TMAX equals CSHIFT and TMIN equals CSHIFT+L1, the length of the first half-period. Here is an R function for illustration: f <- function( camp=10, # circadian amplitude cshift=3, # circadian shift, time of occurrence of maximum time=seq(0, 100, length=100), # time axis L1=5 # length of first interval or time of occurrence of the minimum (time after cshift) ) { L2 <- 24-L1 # derived length of second period # Split time axis from (0, L1) and (L1, L2) to (0, 12) and (12, 24), # then stretch the two intervals individually (note that time gets stretched) t24 <- (time-cshift)%%24 # remainder of division by 24 t1 <- t24/L1*12 # part 1 from (0, L1) to (0, 12) t2 <- (t24-L1)/L2*12 + 12 # part 2 from (L1, L2) to (12, 24) tx <- ifelse(t24 < L1, t1, t2) # step-function but at L1, t1=t2=12 and gradient 0 y <- camp*cos(tx/24*2*pi) # output some of the data and a graph data <- cbind(time, t24, t1, t2, tx, y) print(data[1:min(25, nrow(data)), ]) str <- paste("L1=", L1, ", L2=", L2, ", cshift=", cshift, ", camp=", camp, sep="") plot (time, y, type="l", lwd=5, main=str) abline(v=seq( cshift, max(time), by=24), col="darkblue") # daily max abline(v=seq(L1+cshift, max(time), by=24), col="blue") # daily min } # Call the function: f() # default arguments f(cshift=10, L1=12) # standard with equal lengths f(cshift=10, L1=6) f(cshift=10, L1=6, time=seq(0, 25, length=1000)) # zoom in ----- Andreas Krause, PhD Director, Lead Scientist Modeling and Simulation Dept. of Clinical Pharmacology Actelion Pharmaceuticals Ltd, Hegenheimermattweg 91, CH-4123 Allschwil, Switzerland