-
Notifications
You must be signed in to change notification settings - Fork 0
Open
Description
When trying to fit individual data using els method, I face this problem. Any particular suggestions?
library(dplyr)
library(mrgsolve)
library(nloptr)
library(ggplot2)
library(rlang)
library(klava)
data <- read_csv("data.csv")
mod <- mread("model1.cpp")
data <- data %>% filter(CMT %in% c(1, 2)) %>% mutate(AMT = replace(AMT, AMT==".", 0))
data <- data %>%
mutate(EVID = if_else(CMT==2, 0, 1))
ggplot(data, aes(TIME, DV))+ geom_point()+theme_bw()
theta <- all_ident(tvcl = 18, tvv = 30, tvkatr = 6.5,
tvic50 = 0.22,
tvcort_in = 1.66,
a0 = 0.44,
a1 = 0.22,
b1 = 0.16,
a2 = 0.04,
b2 = -0.05,
pi = 3.14159265359,
hill = 3)
data <- data %>% mutate(AMT = as.numeric(AMT))
data <- data %>% mutate(CMT = if_else(CMT==2, 0, 1))
fit <- fit_nl(theta, data, mod = mod, pred_name= "Conc", cov_step=F,
pred_initial=TRUE, ofv = els)
Error message:
Fitting with els ...Error in res^2/sigma : non-numeric argument to binary operator
the data file contains: ID, TIME, DV, CMT, MDV, WT
for the model file
[init]
gut = 0 ,
cent = 0.44,
endo =0.44,
trans1 = 0,
trans2 = 0,
trans3 =0
[param]
tvcl =22.6
tvv = 41.2
tvkatr = 10.5
a0 = 0.44
a1 = 0.22
b1 = 0.16
a2 = 0.04
b2 = -0.05
pi = 3.14159265359
hill = 3
tvic50 = 1.02
tvcort_in = 1.66
wt = 20
[omega]
0.0818 0.0431 0.15 0.111 0.941
[sigma]
0.0326 0.0098
[main]
double cl = tvcl*pow((wt/70), 0.75);//*exp(ETA(1));
double v = tvv*(wt/70);//*exp(ETA(2));
double kel = cl/v;
double katr = tvkatr;//*exp(ETA(3));
double ic50 = tvic50;//*exp(ETA(4));
double cort_in = tvcort_in*pow((wt/70), 0.75 );//* exp(ETA(5));
double s2 = v/100;
[ode]
#define t SOLVERTIME
double rcort = a0*kel+((a1*kel+b1*(2*pi*1/24))*cos(2*pi*t/24) + (b1*kel+a1*(2*pi/24))*sin(2*pi*t/24) + (a2*kel+b2*(2*pi*1/12))*cos(2*pi*t/12) + (b2*kel+a2*(2*pi/12))*sin(2*pi*t/12));
if(rcort<=0) rcort = 0.001;
double cp = cent/s2 ;
double effect = pow(cp,hill)/(pow(ic50,hill)+pow(cp,hill));
//double effohp = pow(cp,hill)/(pow(i50,hill)+pow(cp,hill));
//double effd4a = pow(cp,hill)/(pow(d50,hill)+pow(cp,hill));
//------------ Differential Equations-----------------//
dxdt_gut = -katr*gut;
dxdt_trans1 = katr*(gut-trans1);
dxdt_trans2 = katr*(trans1-trans2);
dxdt_trans3 = katr*(trans2-trans3);
dxdt_cent = katr*trans3 - kel*cent;
dxdt_endo = cort_in*rcort*(1-effect)-kel*endo;
//dxdt_ohp = rohp*rcort*(1-effohp)-kout_ohp*ohp;
//dxdt_d4a = rd4a*rcort*(1-effd4a)-kout_d4a*d4a;
[table]
double Conc = (cp + endo)* (1+EPS(1)) + EPS(2);
//double OHP = ohp;//*(1+EPS(3));
//if(OHP <= 0) OHP = 0.001;
//double D4a = d4a;//*(1+EPS(4));
//if(D4a <= 0) D4a = 0.001;
if(Conc <= 0) Conc=0.001;
[capture]
Conc
Metadata
Metadata
Assignees
Labels
No labels