Git Product home page Git Product logo

nlmixr's Introduction

nlmixr: an R package for population PKPD modeling

NOTE: nlmixr is curently archived on CRAN. You can also use nlmixr2 which is on CRAN currenty. See http://nlmixr2.org or https://github.com/nlmixr2/nlmixr2

R build status Codecov test coverage CRAN checks CRAN total downloads CRAN monthly downloads


nlmixr

nlmixr is an R package for fitting general dynamic models, pharmacokinetic (PK) models and pharmacokinetic-pharmacodynamic (PKPD) models in particular, with either individual data or population data. The nlme and SAEM estimation routines can be accessed using a universal user interface (UUI), that provides universal model and parameter definition syntax and results in a fit object that can be used as input into the Xpose package. Running nlmixr using the UUI is described in this vignette.

Under the hood nlmixr has five main modules:

  1. dynmodel() and its mcmc cousin dynmodel.mcmc() for nonlinear dynamic models of individual data;
  2. nlme_lin_cmpt()for one to three linear compartment models of population data with first order absorption, or i.v. bolus, or i.v. infusion using the nlme algorithm;
  3. nlme_ode() for general dynamic models defined by ordinary differential equations (ODEs) of population data using the nlme algorithm;
  4. saem_fit for general dynamic models defined by ordinary differential equations (ODEs) of population data by the Stochastic Approximation Expectation-Maximization (SAEM) algorithm;
  5. gnlmm for generalized non-linear mixed-models (possibly defined by ordinary differential equations) of population data by the adaptive Gaussian quadrature algorithm.

A few utilities to facilitate population model building are also included in nlmixr.

Documentation can be found at https://nlmixrdevelopment.github.io/nlmixr/, and we maintain a comprehensive and ever-growing guide to using nlmixr at our bookdown site.

More examples and the associated data files are available at https://github.com/nlmixrdevelopment/nlmixr/tree/master/vignettes.

We recommend you have a look at RxODE, the engine upon which nlmixr depends, as well as xpose.nlmixr, which provides a link to the seminal nonlinear mixed-effects model diagnostics package xpose, and shinyMixR, which provides a means to build a project-centric workflow around nlmixr from the R command line and from a streamlined shiny front-end application. Members of the nlmixr team also contribute to the ggPMX, xgxr and pmxTools packages. For PKPD modeling (with ODE and dosing history) with Stan, check out Yuan Xiong's package PMXStan.

Installation

When on CRAN, you can install the released version of nlmixr from CRAN with:

install.packages("nlmixr")

And the development version from GitHub with:

# install.packages("devtools")
devtools::install_github("nlmixrdevelopment/nlmixr")

nlmixr's People

Contributors

baltcir1 avatar bgoodri avatar billdenney avatar iamstein avatar jranke avatar kestrel99 avatar kiranmaiganji avatar lionel- avatar mattfidler avatar mntrame avatar rikschoemaker avatar vupil avatar wwang-at-github avatar

Stargazers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

Watchers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

nlmixr's Issues

Theophylline PK parameter estimates

Hello,
First, thank you for the terrific work in developing this package!
I’m fitting the Theophylline concentration data using nlmixr and comparing the results to the ones I obtain with saemix or nlme. The parameter estimates are similar between saemix and nlme, but differ (too much to be due to approx. or algorithmic differences) with nlmixr – see the estimates for lCl and lV in the attachment.
I suspect it is a simple variable transformation which is too obvious for me to spot (sic), but I can’t find it. Your insights on this would be much appreciated.
Thanks a lot!

Snippet of nlme code:

startvec1<-c(lKa=0.5, lCl=0.75, lV=3.45)
nform<- ~(Dose*exp(lKa)*(exp(-(exp(lCl)/exp(lV))*Time)-exp(-exp(lKa)*Time)))/(exp(lV)*(exp(lKa)-(exp(lCl)/exp(lV))))
nlme.theomod<-deriv(nform, namevec=c("lKa", "lCl", "lV"), function.arg=c("Dose", "Time", "lKa", "lCl", "lV"))
Theo.nlme<-nlme(Concentration~nlme.theomod(Dose, Time, lKa, lCl, lV), 
                data=groupedData(Concentration~Time|Id, data=theodf), 
                fixed=list(lKa~1, lCl~1, lV~1), random=pdDiag(lKa+lCl+lV~1), start=startvec1)

Snippets of nlmixr code:

uif1 <- function() {
  ini({
    lka <- 0.45
    lCl <- 1
    lV  <- 3.45
    eta.lka ~ 0.6
    eta.lCl ~ 0.3
    eta.lV  ~ 0.1
    add.err <- 0.7
  })
  model({
    ka <- exp(lka+eta.lka)
    cl <- exp(lCl+eta.lCl)
    v <- exp(lV+eta.lV)
    linCmt() ~ add(add.err)
  })
}
nlmxir.fit1<-nlmixr(uif1, indf, est="nlme", calc.resid=FALSE)

TheoCompare.pdf

dll needs to be removed with nlmxr/SAEM before a new model can be run

When you create an SAEM model, e.g. using the following syntax:

saem_fit <- gen_saem_user_fn(model=lincmt(ncmt=1, oral=FALSE,infusion=TRUE))

and then want to run a new model, you need to explicitly remove the dll before this works using:

dyn.unload("saem_main.dll")
saem_fit <- gen_saem_user_fn(model=lincmt(ncmt=2, oral=FALSE,infusion=TRUE))

Although it is nice that the dll is persistent so it can be reused for new data with the same model, it would be much nicer if the dll can be referenced explicitly (like with nlme) and put in a separate subdirectory instead of being created in the run directory. This may also be related to the issue with running multiple SAEM instances in parallel (see other issue).

nlmeControl documentation issues

I am interested in seeing the likelihood (or log-likelihood) for a combination of model, known parameter values, and a data set; this would be the value prior to any parameter optimization. I've read the nlmixr documentation fairly carefully but it's not yet clear to me how best to do this. I assume that this should be controlled by an option in nlmeControl or saemControl (but only the latter has an entry in the documentation [let's call this sub-issue 1]). In the Package 'nlmixr' documentation (2017-11-09) under gnlmm, in the Examples section, there is this line:

fit = gnlmm(llik, pump, inits,
control=list(
reltol.outer=1e-4,
optim.outer="nmsimplex",
nAQD=5
)
)

What does nAQD control? [sub-issue 2] It would be nice to have some discussion or reference to the control parameters.

Would something like this be suitable for getting only the initial objective function value or log-likelihood?

fit <- nlmixr(one.compartment.model, data, est="nlme",
              control = nlmeControl(maxIter = 1, pnlsMaxIter = 1, msMaxIter = 1))

Is this overkill or the best way to get this information? [sub-issue 3]

Thank-you,
James

source code for Acop poster

do you have a github repository for the simulation and estimation models for the ACoP poster in 2016? Curious about results of the sparse datasets, so it would be great if you could provide a repo of all the work that generated the datasets and the model files for nonmem and nlmixr

VPC does not work for NLME objects

The vpc() function does not work for fit objects generated by nlme_lin_cmt(). Fails with error:

Error in plot.window(...) : need finite 'ylim' values

Error installing on Windows 10

Hello,

I tried to following the installing instruction here:
https://github.com/nlmixrdevelopment/nlmixr/blob/master/inst/Installing_nlmixr.pdf

But as soon as I get to the step:
install_github('nlmixrdevelopment/nlmixr')
Or:
install_github('nlmixrdevelopment/RxODE')

I keep seeing this error message:

  • installing source package 'nlmixr' ...
    Warning in file(file, if (append) "a" else "w") :
    cannot open file 'C:/Users/Name/Documents/R/win-library/3.4/nlmixr/DESCRIPTION': No such file or directory
    Error in file(file, if (append) "a" else "w") :
    cannot open the connection
    ERROR: installing package DESCRIPTION failed for package 'nlmixr'
  • removing 'C:/Users/Name/Documents/R/win-library/3.4/nlmixr'
    Installation failed: Command failed (1)

I tried running R in administrator mode but that didn't help. Does anyone know how I can solve this?

Thank you.

AMT in theo_sd and theo_md are not the doses administered.

I used the theophylline data as given by theo_sd from nlmixr and fit them in NONMEM. (Unfortunately I still need to use NONMEM because of its MAXEVALS= 0 feature; when Issue #32 is resolved, please let me know!) I took the AMT column at face value and got a THETA vector, no problem. I then put this THETA vector into Madonna but used the known amount of the dose for my Madonna code, 320 mg. The curve was way off from the data. This was puzzling until I realized that the AMT column in the theo_sd data is wrong; the actual amount should be AMT * WT. The AMT is really normalized by weight: if you take AMT times WT you get very close to 320 (mg) for all patients except one. In the original experiment 12 volunteers were given a dose of 320 mg. I don't know why patient 9 is off from that; his value of AMT*WT is 267.840 mg.

If you look at the documentation of these data just from R, ?Theoph says that the Dose column is normalized for weight, but when the Theoph dataset was made more NONMEM-like by the addition of EVID and CMT and renamed columns, apparently Dose was taken as AMT directly.

Thanks for your efforts!

issue with installation in window

When running install_github("nlmixrdevelopment/nlmixr"), ar command is not found, wondering what's wrong with that. Below is the error message.

ar crv libdparse.a arg.o parse.o scan.o symtab.o util.o dparse_tree.o
make[1]: ar: Command not found
make[1]: *** [libdparse.a] Error 127
make[1]: Leaving directory `/cygdrive/c/Users/Dark/AppData/Local/Temp/RtmpC8oA2x/devtoolsf4c33634b58/nlmixrdevelopment-nlmixr-45fd1f4/src/d'
make: *** [libdparse] Error 2
Warning: running command 'make -f "Makevars" -f "C:/PROGRA1/R/R-341.1/etc/x64/Makeconf" -f "C:/PROGRA1/R/R-341.1/share/make/winshlib.mk" SHLIB_LDFLAGS='$(SHLIB_CXXLDFLAGS)' SHLIB_LD='$(SHLIB_CXXLD)' SHLIB="nlmixr.dll" WIN=64 TCLBIN=64 ' had status 2
ERROR: compilation failed for package 'nlmixr'

  • removing 'C:/Users/Dark/Documents/R/win-library/3.4/nlmixr'
    Installation failed: Command failed (1)

want to use nlmixr on bootstrapped data

Below is my session output. I made a function, funBS(), that generates bootstrapped data from a data frame that has ID as its first column. However, when I try to use the bootstrapped data with nlmixr, I get an error and no object is created. Note that much of this is output from what I just pasted in for Issue #28 but here is the full info. The bootstrap part begins with funBS(). Thank-you for your help.

R version 3.4.3 (2017-11-30) -- "Kite-Eating Tree"
Copyright (C) 2017 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> # Preliminaries
> setwd("C:/Users/James/Documents/pharmacometrics/testing")
> library(nlmixr)
Loading required package: nlme
Loading required package: RxODE
> # rxVersion(echo=TRUE)
> rxClean()
[1] TRUE
> rxVersion(echo=TRUE)
  _____         ____  _____  ______
 |  __ \       / __ \|  __ \|  ____| 0.7-0
 | |__) |__  _| |  | | |  | | |__
 |  _  / \ \/ / |  | | |  | |  __|
 | | \ \  >  <| |__| | |__| | |____
 |_|  \_\/_/\_\\____/|_____/|______|
> # Use some pre-loaded data.
> str(theo_sd)  # like Theoph but with ID instead of Subject, also other NONMEM-like
'data.frame':	144 obs. of  7 variables:
 $ ID  : int  1 1 1 1 1 1 1 1 1 1 ...
 $ TIME: num  0 0 0.25 0.57 1.12 2.02 3.82 5.1 7.03 9.05 ...
 $ DV  : num  0 0.74 2.84 6.57 10.5 9.66 8.58 8.36 7.47 6.89 ...
 $ AMT : num  4.02 0 0 0 0 0 0 0 0 0 ...
 $ EVID: int  101 0 0 0 0 0 0 0 0 0 ...
 $ CMT : int  1 2 2 2 2 2 2 2 2 2 ...
 $ WT  : num  79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 ...
> # 'uif' for Unified Interface Function, 1 compartment, lin for linCmt()
> uif.1cm.lin <- function() {
+   ini({
+     logtcl <- -3.2      # log typical value Cl (L/hr)
+     logtka <-  0.5      # log typical value Ka (1/hr)
+     logtv  <- -1        # log typical value V  (L)
+     
+     # error model
+     add.err <- 0.1
+     
+     # Initial estimates for IIV variances
+     # Labels work for single parameters.
+     eta.cl ~ 2
+     eta.ka ~ 1    ## BSV Ka
+     eta.v  ~ 1
+   })
+   model({
+     cl <- exp(logtcl + eta.cl)
+     ka <- exp(logtka + eta.ka)
+     v  <- exp(logtv  + eta.v)
+     linCmt() ~ add(add.err)
+   })
+ }
> # NLME fittings
> 
> fit.nlme.1cm.lin.theo_sd     <- nlmixr(uif.1cm.lin, theo_sd, est="nlme")

**Iteration 1
LME step: Loglik: -182.2318, nlminb iterations: 1
reStruct  parameters:
      ID1       ID2       ID3 
1.0022617 0.2783009 2.0758984 
PNLS step: RSS =  63.14922 
 fixed effects: -3.211937  0.4479979  -0.7859318  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
   fixed reStruct 
0.272375 3.267309 

**Iteration 2
LME step: Loglik: -179.291, nlminb iterations: 9
reStruct  parameters:
       ID1        ID2        ID3 
0.96385527 0.08333435 1.63552295 
PNLS step: RSS =  63.28221 
 fixed effects: -3.211535  0.4441505  -0.7863391  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.008662381 0.172135330 

**Iteration 3
LME step: Loglik: -179.3381, nlminb iterations: 8
reStruct  parameters:
       ID1        ID2        ID3 
0.96129928 0.07085834 1.63885629 
PNLS step: RSS =  63.2196 
 fixed effects: -3.211732  0.4458486  -0.7861743  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.003808629 0.082919683 

**Iteration 4
LME step: Loglik: -179.3201, nlminb iterations: 7
reStruct  parameters:
       ID1        ID2        ID3 
0.96243507 0.07618717 1.63736945 
PNLS step: RSS =  63.22991 
 fixed effects: -3.211823  0.4451494  -0.7862432  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.001570583 0.013884126 

**Iteration 5
LME step: Loglik: -179.3277, nlminb iterations: 4
reStruct  parameters:
       ID1        ID2        ID3 
0.96197111 0.07396744 1.63796900 
PNLS step: RSS =  63.24112 
 fixed effects: -3.211722  0.4454457  -0.7862145  
 iterations: 6 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.0006651255 0.0147544693 

**Iteration 6
LME step: Loglik: -179.3246, nlminb iterations: 1
reStruct  parameters:
       ID1        ID2        ID3 
0.96215433 0.07487426 1.63772430 
PNLS step: RSS =  63.24112 
 fixed effects: -3.211722  0.4454457  -0.7862145  
 iterations: 1 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.000000e+00 3.910435e-10 
Using sympy via SnakeCharmR
## Calculate ETA-based prediction and error derivatives:
Calculate Jacobian...................done.
Calculate sensitivities.......
done.
## Calculate d(f)/d(eta) 
## ...
## done
## ...
## done
C:/RTools/3.4/mingw_64/bin/gcc  -I"C:/R/R-34~1.3/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_ebf0e5096f80d0203de27cbc48be0f26_x64.c -o rx_ebf0e5096f80d0203de27cbc48be0f26_x64.o
C:/RTools/3.4/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_ebf0e5096f80d0203de27cbc48be0f26_x64.dll tmp.def rx_ebf0e5096f80d0203de27cbc48be0f26_x64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
C:/RTools/3.4/mingw_64/bin/gcc  -I"C:/R/R-34~1.3/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_40eae078bd4f9e0651450c3b310d5e9a_x64.c -o rx_40eae078bd4f9e0651450c3b310d5e9a_x64.o
C:/RTools/3.4/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_40eae078bd4f9e0651450c3b310d5e9a_x64.dll tmp.def rx_40eae078bd4f9e0651450c3b310d5e9a_x64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
The model-based sensitivities have been calculated.
It will be cached for future runs.
Diagonal form: sqrt
     [,1]
[1,]    1
Calculate symbolic inverse...done
Calculate symbolic determinant of inverse...done
Calculate d(Omega)/d(Est) and d(Omega^-1)/d(Est)...
.0
...done.
Calculating Table Variables...
done
Warning message:
In nlmixrUI.nlme.var(obj) :
  Initial condition for additive error ignored with nlme
nlmixr UI combined dataset and properties
 $ par.hist         : Parameter history (if available)
 $ par.hist.stacked : Parameter history in stacked form for easy plotting (if available)
 $ par.fixed        : Fixed Effect Parameter Table
 $ eta              : Individual Parameter Estimates
 $ seed             : Seed (if applicable)
 $ model.name       : Model name (from R function)
 $ data.name        : Name of R data input
> # Only 1 iteration for DFN.
> fit.nlme.1cm.lin.theo_sd.init <- nlmixr(uif.1cm.lin, theo_sd, est="nlme",
+                                         control = nlmeControl(maxIter = 1,
+                                                               pnlsMaxIter = 1,
+                                                               msMaxIter = 1)
+                                         )

**Iteration 1
LME step: Loglik: -182.2318, nlminb iterations: 1
reStruct  parameters:
      ID1       ID2       ID3 
1.0022617 0.2783009 2.0758984 
PNLS step: RSS =  87.52968 
 fixed effects: -3.230002  0.5346704  -0.7523971  
 iterations: 1 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
3.290853e-01 4.201928e-06 
Error in nlme.formula(model = DV ~ (nlmixr::nlmeModList("user_fn"))(logtcl,  : 
  maximum number of iterations (maxIter = 1) reached without convergence
In addition: Warning message:
In nlmixrUI.nlme.var(obj) :
  Initial condition for additive error ignored with nlme
> # No object created.
> # So try to expand the tolerance:
> fit.nlme.lin.theo_sd.init <- nlmixr(uif.1cm.lin, theo_sd, est="nlme",
+                                     control = nlmeControl(maxIter = 1,
+                                                           pnlsMaxIter = 1,
+                                                           msMaxIter = 1),
+                                     tolerance = 1 # deliberately a huge tolerance
+                                     )
Error in (function (model, data = sys.frame(sys.parent()), fixed, random = fixed,  : 
  unused argument (tolerance = 1)
In addition: Warning message:
In nlmixrUI.nlme.var(obj) :
  Initial condition for additive error ignored with nlme
> # No object created.
> # Just curious, will nlmixrUI(fun) work? It's been removed from latest documentation.
> nlmixrUI(uif.1cm.lin)
## 1-compartment model with first-order absorption in terms of Cl
## Initialization:
################################################################################
Fixed Effects ($theta):
logtcl logtka  logtv 
  -3.2    0.5   -1.0 

Omega ($omega):
     [,1] [,2] [,3]
[1,]    2    0    0
[2,]    0    1    0
[3,]    0    0    1

## Model:
################################################################################
    cl <- exp(logtcl + eta.cl)
    ka <- exp(logtka + eta.ka)
    v  <- exp(logtv  + eta.v)
> # Seems to still work, not sure why its documentation was removed from nlmixr.pdf.
> 
> # Can nlmixr work on bootstrapped data?
> funBS <- function(df, B) {
+   # df is a data frame with column 1 called ID.
+   # B is the number of bootstrapped data sets.
+   subjects <- sort(as.integer(unique(df$ID)))
+   
+   # Initialize an empty data frame to hold bootstrapped values.
+   BsIteration.and.Subject <- data.frame(BsIteration=numeric(),
+                                         ID=character(),
+                                         stringsAsFactors=FALSE)
+   temp <- setNames(data.frame(matrix(ncol = 2, nrow = 0), stringsAsFactors=FALSE),
+                    c("BsIteration", "ID"))
+   for (b in 1:B){
+     # Create a vector of bootstrapped subject IDs.
+     BSsubjects <- sample(1:length(subjects), length(subjects), replace=TRUE)
+     for (i in BSsubjects){
+       # The next line relies on ID being column 1. Index b gets recycled.
+       temp <- as.data.frame( c(b, df[df$ID==i,][1]), stringsAsFactors=FALSE)
+       names(temp) <- c("BsIteration", "ID")
+       BsIteration.and.Subject <- rbind(BsIteration.and.Subject, temp)
+     }
+   }
+   # Create the rest of the empty vector (to avoid lots of rbind on lots of data).
+   other <- setNames(data.frame(matrix(nrow=nrow(BsIteration.and.Subject),
+                                       ncol=ncol(df[-1])),
+                                stringsAsFactors=FALSE),
+                     names(df[-1]) )
+   # 'BSDS' for bootstrapped data set
+   BSDS <- cbind(BsIteration.and.Subject, other)
+   
+   # Fill in BSDS.
+   for (i in BSDS$ID){
+     BSDS[,c(-1,-2)] <- df[df$ID==i,2:ncol(df)]
+   }
+   
+   return(list(BSsubjects, BSDS))
+ }
> BSDS002 <- funBS(theo_sd, 2)  # 002 suffix is number of BS samples
> BSDS002[[1]]
 [1] 10 12  8 11  1  7  5 10  5  9  2  7
> BSDS002[[2]]$ID
  [1]  2  2  2  2  2  2  2  2  2  2  2  2 11 11 11 11 11 11 11 11 11 11 11 11 10 10 10
 [28] 10 10 10 10 10 10 10 10 10  2  2  2  2  2  2  2  2  2  2  2  2  9  9  9  9  9  9
 [55]  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  6  6  6  6  6  6  6  6  6
 [82]  6  6  6  9  9  9  9  9  9  9  9  9  9  9  9  1  1  1  1  1  1  1  1  1  1  1  1
[109] 10 10 10 10 10 10 10 10 10 10 10 10  4  4  4  4  4  4  4  4  4  4  4  4  7  7  7
[136]  7  7  7  7  7  7  7  7  7 10 10 10 10 10 10 10 10 10 10 10 10 12 12 12 12 12 12
[163] 12 12 12 12 12 12  8  8  8  8  8  8  8  8  8  8  8  8 11 11 11 11 11 11 11 11 11
[190] 11 11 11  1  1  1  1  1  1  1  1  1  1  1  1  7  7  7  7  7  7  7  7  7  7  7  7
[217]  5  5  5  5  5  5  5  5  5  5  5  5 10 10 10 10 10 10 10 10 10 10 10 10  5  5  5
[244]  5  5  5  5  5  5  5  5  5  9  9  9  9  9  9  9  9  9  9  9  9  2  2  2  2  2  2
[271]  2  2  2  2  2  2  7  7  7  7  7  7  7  7  7  7  7  7
> head(BSDS002[[2]]$ID, 15)
 [1]  2  2  2  2  2  2  2  2  2  2  2  2 11 11 11
> head(BSDS002[[2]], 15)
   BsIteration ID  TIME   DV  AMT EVID CMT   WT
1            1  2  0.00 0.00 4.95  101   1 64.6
2            1  2  0.00 0.15 0.00    0   2 64.6
3            1  2  0.25 0.85 0.00    0   2 64.6
4            1  2  0.50 2.35 0.00    0   2 64.6
5            1  2  1.02 5.02 0.00    0   2 64.6
6            1  2  2.02 6.58 0.00    0   2 64.6
7            1  2  3.48 7.09 0.00    0   2 64.6
8            1  2  5.00 6.66 0.00    0   2 64.6
9            1  2  6.98 5.25 0.00    0   2 64.6
10           1  2  9.00 4.39 0.00    0   2 64.6
11           1  2 12.05 3.53 0.00    0   2 64.6
12           1  2 24.22 1.15 0.00    0   2 64.6
13           1 11  0.00 0.00 4.95  101   1 64.6
14           1 11  0.00 0.15 0.00    0   2 64.6
15           1 11  0.25 0.85 0.00    0   2 64.6
> tail(BSDS002[[2]], 15)
    BsIteration ID  TIME   DV  AMT EVID CMT   WT
274           2  2  9.00 4.39 0.00    0   2 64.6
275           2  2 12.05 3.53 0.00    0   2 64.6
276           2  2 24.22 1.15 0.00    0   2 64.6
277           2  7  0.00 0.00 4.95  101   1 64.6
278           2  7  0.00 0.15 0.00    0   2 64.6
279           2  7  0.25 0.85 0.00    0   2 64.6
280           2  7  0.50 2.35 0.00    0   2 64.6
281           2  7  1.02 5.02 0.00    0   2 64.6
282           2  7  2.02 6.58 0.00    0   2 64.6
283           2  7  3.48 7.09 0.00    0   2 64.6
284           2  7  5.00 6.66 0.00    0   2 64.6
285           2  7  6.98 5.25 0.00    0   2 64.6
286           2  7  9.00 4.39 0.00    0   2 64.6
287           2  7 12.05 3.53 0.00    0   2 64.6
288           2  7 24.22 1.15 0.00    0   2 64.6
> # Bootstrapped data set BSDS002 seems reasonable.
> # nlmixr on bootstrapped data unfortunately fails:
> fit.nlme.1cm.lin.BSDS002.full <- nlmixr(uif.1cm.lin, BSDS002, est="nlme")
Error in apply(data, 1, function(x) { : 
  dim(X) must have a positive length
> 

installation problem

installation_tests.txt
failure at step 7.txt

I was following along in "Installation on Windows" (https://github.com/nlmixrdevelopment/nlmixr) and things were basically working okay, although there were a few unit tests that failed (after a long time) in Step 6 (Install RxODE), last sub-step. Then at Step 7, to actually install nlmixr, the other packages got installed fine, but nlmixr compilation failed. I do not know why.

I pasted the output for Steps 1-6 in its own file (which is relatively big because of all the testing in Step 6), and the problematic Step 7 in its own file, both attached.

By the way, as for the Windows installer, it installs its own version of R and was invisible to the version of R I already had running just fine seen by RStudio, so I uninstalled all that prior to following the steps in "Installation on Windows".

Compilation Terminated - RcppArmadillo.h

Hi,

I am trying to install nlmixr but I am getting the following error:

g++ -m64 -std=gnu++0x -I/usr/include/R -DNDEBUG  -I"/usr/lib64/R/library/dparser/include" -I"/usr/lib64/R/library/Rcpp/include" -I"/usr/lib64/R/library/RcppEigen/include" -I"/usr/lib64/R/library/StanHeaders/include" -I"/usr/lib64/R/library/BH/include" -I/usr/local/include  -Id -I../inst/include -DBOOST_DISABLE_ASSERTS -I"`"/usr/lib64/R/bin/Rscript" -e 'cat(file.path(system.file("", package = "RcppArmadillo"),"include"))'`" -fpic  -O2 -g -pipe -Wall -Wp,-D_FORTIFY_SOURCE=2 -fexceptions -fstack-protector --param=ssp-buffer-size=4 -m64 -mtune=generic -c RcppExports.cpp -o RcppExports.o
Rscript execution error: Exec format error
g++ -m64 -std=gnu++0x -I/usr/include/R -DNDEBUG  -I"/usr/lib64/R/library/dparser/include" -I"/usr/lib64/R/library/Rcpp/include" -I"/usr/lib64/R/library/RcppEigen/include" -I"/usr/lib64/R/library/StanHeaders/include" -I"/usr/lib64/R/library/BH/include" -I/usr/local/include  -Id -I../inst/include -DBOOST_DISABLE_ASSERTS -I"`"/usr/lib64/R/bin/Rscript" -e 'cat(file.path(system.file("", package = "RcppArmadillo"),"include"))'`" -fpic  -O2 -g -pipe -Wall -Wp,-D_FORTIFY_SOURCE=2 -fexceptions -fstack-protector --param=ssp-buffer-size=4 -m64 -mtune=generic -c focei_fit.cpp -o focei_fit.o
Rscript execution error: Exec format error
focei_fit.cpp:1:27: fatal error: RcppArmadillo.h: No such file or directory
compilation terminated.
make: *** [focei_fit.o] Error 1
ERROR: compilation failed for package ‘nlmixr’
* removing ‘/usr/lib64/R/library/nlmixr’
Installation failed: Command failed (1)
> q()

I located RcppArmadillo.h and I do have it in my library:

[root@amrvlp000005185 Python-3.6.5]# locate -i RcppArmadillo.h
/usr/lib64/R/library/RcppArmadillo/include/RcppArmadillo.h
/var/lib/rstudio-connect/packrat/3.4.3/v2/library/RcppArmadillo/06d74d305dc7f3c53054e88efd49b097/RcppArmadillo/include/RcppArmadillo.h
/var/lib/rstudio-connect/packrat/3.4.3/v2/library/RcppArmadillo/dd3c11bed5cee97ee28e5543659ed01a/RcppArmadillo/include/RcppArmadillo.h
/var/lib/rstudio-connect/packrat/3.4.3/v2/library/RcppArmadillo/f280a22174c410654514ef1860a246be/RcppArmadillo/include/RcppArmadillo.h

Please advise. Thank you in advance.

allow "+

I sometimes use the following syntax for clarity purposes:

d/dt(A2) = +(Q/V)*A1 - (Q/V2)*A2

the initial (unnecessary) + makes the compilation fail, see error below. Would be nice if RxODE would allow this, or else give back a warning message that this is not allowed.

Compiling RxODE differential equations...Error in dyn.load(finalDll, local = FALSE) : 
  unable to load shared object '<...._.so':
  /home/ron...._.so: cannot open shared object file: No such file or directory
In addition: Warning message:
In nlmixrData.default(data) :
  NONMEM-style data converted to nlmixr/RxODE-style data.
Error in dyn.load(basename(finalDll), local = FALSE) : 
  unable to load shared object '/tmp/RtmppRlnoA/file440473e62dc/rx_fabd2d70cc2a164953c4931d0c425bb3_.so':
  /tmp/RtmppRlnoA/file440473e62dc/rx_fabd2d70cc2a164953c4931d0c425bb3_.so: cannot open shared object file: No such file or directory
Error in rxCompile.character(model, extraC = extraC, debug = debug, calcJac = calcJac,  : 
  Error loading model.

Freezing some of the parameters in SAEM in nlmixr

Hello Rik, Wenping, and all,

I am trying out saem in nlmixr and wondering about how to freeze some of the parameters. I am hoping to read in some parameters directly from the dataset, but it does not work. It runs, but kind of defaults to some theta in estimates (I know for sure that it is not correct estimate at all). Here is what I am trying. I have V1, V2, CL, and Q defined in the dataset (Sims3) (same estimates for all IDs).

ODE <-
"C1 = centr/V1;
C2 = peri/V2;
d/dt(depot) = -Kadepot;
d/dt(centr) = F1
Kadepot - CLC1 - QC1 + QC2;
d/dt(peri) = QC1 - QC2;
"
mypar <- function(lKa, lKdeg, V1, V2, CL, Q)
{Ka = exp(lKa)
F1=F1
V1=V1
V2=V2
CL=CL
Q=Q
}

Sims3$COV <- 1
m2 = RxODE(ODE, modName="m2")
PRED = function() C1
saem_fit <- gen_saem_user_fn(model=m2, PKpars=mypar, pred=PRED)
model = list(saem_mod=saem_fit, res.mod=3, covars="COV")
inits = list(theta=c(Ka=0.1))
cfg = configsaem(model, data=Sims3, inits)
cfg$print = 50
fit = saem_fit(cfg)
summary(fit)
plot(fit)

Instead if I switch mypar with something like following, it works and gives me reliable estimates.
mypar <- function(lKa, V1, V2, CL, Q)
{Ka = exp(lKa)
V1=2
V2=3
CL=0.3
Q=0.2
F1=0.65
}

Any way to define parameters that I want to estimate in dataset rather than in function?

Thanks,
Krina

error compiling for Windows - solved

R 3.3.2, Rtools 34, Windows 8.1 (x64):

d/libdparse.a: error adding symbols: Archive has no index; run ranlib to add one
collect2.exe: error: ld returned 1 exit status
no DLL was created
ERROR: compilation failed for package 'nlmixr'

SAEM syntax requires a defined covariate

Currently, SAEM syntax will only run when a covariate is explicitly specified:

The work-around is to specify a dummy covariate but that should not be necessary:

saem_fit <- gen_saem_user_fn(model=lincmt(ncmt=1, oral=FALSE))
nmdat$WT<-1
model = list(saem_mod=saem_fit, res.mod=2,covars="WT")
etc etc

nlme_ode running error

Hi, thank you for developing such a great package! I am learning your package for my modeling project. I've installed the Github version of nlmixr and RxODE. Could you please look at the error message below and tell me what was wrong? It's nlme_ode() using the theo_md example. The dynmodel() works fine, so I guess it is not an installation problem. Although I did get an error message about the xpose.nlme extension, I read from another post that this is not a problem. Another question is: how can we do naive pool fitting of multiple-dose data (each individual has different dose)? Thank you!

> library(nlmixr)
> ode <- "
+ d/dt(depot) =-KA*depot;
+ d/dt(centr) = KA*depot - KE*centr;
+ "
> dat <- read.table(system.file("examples/theo_md.txt", package = "nlmixr"), head=TRUE)
> mypar <- function(lKA, lKE, lCL)
+ {
+     KA=exp(lKA)
+     KE=exp(lKE)
+     CL=exp(lCL)
+     V = CL/KE
+ }
> 
> specs <- list(fixed=lKA+lKE+lCL~1, random = pdDiag(lKA+lCL~1),
+               start=c(lKA=0.5, lKE=-2.5, lCL=-3.2))
> fit <- nlme_ode(dat, model=ode, par_model=specs, par_trans=mypar,
+                 response="centr", response.scaler="V")
C:/Rtools/mingw_64/bin/gcc  -I"C:/PROGRA~1/R/R-35~1.0/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_97ee19c2b05e8b4a202d9dbdfed6d81d_x64.c -o rx_97ee19c2b05e8b4a202d9dbdfed6d81d_x64.o
C:/Rtools/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_97ee19c2b05e8b4a202d9dbdfed6d81d_x64.dll tmp.def rx_97ee19c2b05e8b4a202d9dbdfed6d81d_x64.o -LC:/PROGRA~1/R/R-35~1.0/bin/x64 -lRblas -LC:/PROGRA~1/R/R-35~1.0/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/PROGRA~1/R/R-35~1.0/bin/x64 -lR
Error in nlme.formula(model = DV ~ (nlmixr::nlmeModList("user_fn"))(lKA,  : 
  step halving factor reduced below minimum in PNLS step
> 

saem + full block / nlme

Hi all,

Great work on the recent developments and release! Got it to work on all platforms now, nice job on that. Basic examples run fine, but encountered a few things when exploring one of the examples a bit deeper, not sure if I'm doing something wrong... Using Oral_1CPT data, and model example as given in help:

  1. when using a full omega matrix (and using SAEM), a crash occurs:

Omega matrix used:

    eta_CL + eta_V + eta_KA ~ c(0.1,
                                0.05, 0.1,
                                0.01, 0.01, 0.1)

Error:

filed95537d84b6b.cpp:21:1001: error: no matching function for call to 'R_pow_di'
/Library/Frameworks/R.framework/Resources/include/Rmath.h:375:8: note: candidate function not viable: requires 2 arguments, but 1 was provided
double R_pow_di(double, int);
       ^
  1. using the same oral 1cmt model with diagonal matrix but estimating with nlme, I get a crash after 4th iteration:

call:

fit   <- nlmixr(object = f, 
                data = d, 
                est = "nlme")

error:

**Iteration 3
LME step: Loglik: -37499.72, nlminb iterations: 23
reStruct  parameters:
       ID1        ID2        ID3 
-0.2728802 -0.3364310 -0.5006904 
varStruct  parameters:
    const 
-2.609917 
Error in nlme.formula(model = DV ~ (nlmeModList("user_fn"))(lCL, eta_CL,  : 
  step halving factor reduced below minimum in PNLS step

Any ideas on these issues? First one seems some matrix definiteness issue perhaps? But second one occurs also with diagonal block.

Thanks!

very different SAEM results and other tutorial issues

The tutorial (https://github.com/nlmixrdevelopment/nlmixr/blob/master/vignettes/running_nlmixr.Rmd) has code for compiling a model for closed form solution using linCmt() and for ODE solution, and also shows how one might use NMLE and SAEM algorithms. With the single dose theophylline data at hand, this gives 4 possibilities. I followed along the tutorial but made a few changes for clarity and for comparison and found that the NMLE result for closed form and for ODE agree very well (but not exactly), whereas the SAEM objective function values, AIC, BIC, and log-likelihood are very different for the explicit versus the ODE models and do not agree with the NMLE results either. The SAEM results seem to be way off. Below is the code and results (and by the way, is there a good way to make all of the following pasted text appear in a uniform typewriter font?):
tutorial.txt
tutorial output.txt

`The tutorial (https://github.com/nlmixrdevelopment/nlmixr/blob/master/vignettes/running_nlmixr.Rmd) has code for compiling a model for closed form solution using linCmt() and for ODE solution, and also shows how one might use NMLE and SAEM algorithms. With the single dose theophylline data at hand, this gives 4 possibilities. I followed along the tutorial but made a few changes for clarity and for comparison and found that the NMLE result for closed form and for ODE agree very well (but not exactly), whereas the SAEM objective function values, AIC, BIC, and log-likelihood are very different for the explicit versus the ODE models and do not agree with the NMLE results either. The SAEM results seem to be way off. Below is the code and results:

R version 3.4.3 (2017-11-30) -- "Kite-Eating Tree"
Copyright (C) 2017 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> setwd("C:/Users/James/Documents/pharmacometrics/practice/nlmixr")
> #https://github.com/nlmixrdevelopment/nlmixr/blob/master/vignettes/running_nlmixr.Rmd
> 
> # Preliminaries
> library(ggplot2)
> library(nlmixr)
Loading required package: nlme
Loading required package: RxODE
> str(theo_sd)
'data.frame':	144 obs. of  7 variables:
 $ ID  : int  1 1 1 1 1 1 1 1 1 1 ...
 $ TIME: num  0 0 0.25 0.57 1.12 2.02 3.82 5.1 7.03 9.05 ...
 $ DV  : num  0 0.74 2.84 6.57 10.5 9.66 8.58 8.36 7.47 6.89 ...
 $ AMT : num  4.02 0 0 0 0 0 0 0 0 0 ...
 $ EVID: int  101 0 0 0 0 0 0 0 0 0 ...
 $ CMT : int  1 2 2 2 2 2 2 2 2 2 ...
 $ WT  : num  79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 79.6 ...
> ggplot(theo_sd, aes(TIME, DV)) +
+   geom_line(aes(group=ID), col="red") + 
+   scale_x_continuous("Time (h)") +
+   scale_y_continuous("Concentration") +
+   labs(title="Theophylline single-dose",
+        subtitle="Concentration vs. time by individual")
> # Try fitting a simple one-compartment PK model to this small dataset.
> # 'uif' for Unified Interface Function, lin for linCmt()
> uif.lin <- function() {
+   ini({
+     logtcl <- -3.2      # log typical value Cl (L/hr)
+     logtka <-  0.5      # log typical value Ka (1/hr)
+     logtv  <- -1        # log typical value V  (L)
+     
+     # error model
+     add.err <- 0.1
+     
+     # Initial estimates for IIV variances
+     # Labels work for single parameters.
+     eta.cl ~ 2
+     eta.ka ~ 1    ## BSV Ka
+     eta.v  ~ 1
+   })
+   model({
+     cl <- exp(logtcl + eta.cl)
+     ka <- exp(logtka + eta.ka)
+     v  <- exp(logtv  + eta.v)
+     linCmt() ~ add(add.err)
+   })
+ }
> # We can alternatively express the same model by ODEs:
> uif.ode <- function() {
+   ini({
+     logtcl <- -3.2      # log typical value Cl (L/hr)
+     logtka <-  0.5      # log typical value Ka (1/hr)
+     logtv  <- -1        # log typical value V  (L)
+     
+     # error model
+     add.err <- 0.1
+     
+     # Initial estimates for IIV variances
+     # Labels work for single parameters.
+     eta.cl ~ 2
+     eta.ka ~ 1    ## BSV Ka
+     eta.v  ~ 1
+   })
+   model({
+     cl <- exp(logtcl + eta.cl)
+     ka <- exp(logtka + eta.ka)
+     v  <- exp(logtv  + eta.v)
+     d/dt(depot)  = -ka * depot
+     d/dt(center) =  ka * depot - cl / v * center
+     cp = center / v
+     cp ~ add(add.err)
+   })
+ }
> # NLME fittings
> 
> fit.nlme.lin <- nlmixr(uif.lin, theo_sd, est="nlme", calc.resid=FALSE)

**Iteration 1
LME step: Loglik: -182.2318, nlminb iterations: 1
reStruct  parameters:
      ID1       ID2       ID3 
1.0022617 0.2783009 2.0758984 
PNLS step: RSS =  63.14922 
 fixed effects: -3.211937  0.4479979  -0.7859318  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
   fixed reStruct 
0.272375 3.267309 

**Iteration 2
LME step: Loglik: -179.291, nlminb iterations: 9
reStruct  parameters:
       ID1        ID2        ID3 
0.96385527 0.08333435 1.63552295 
PNLS step: RSS =  63.28221 
 fixed effects: -3.211535  0.4441505  -0.7863391  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.008662381 0.172135330 

**Iteration 3
LME step: Loglik: -179.3381, nlminb iterations: 8
reStruct  parameters:
       ID1        ID2        ID3 
0.96129928 0.07085834 1.63885629 
PNLS step: RSS =  63.2196 
 fixed effects: -3.211732  0.4458486  -0.7861743  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.003808629 0.082919683 

**Iteration 4
LME step: Loglik: -179.3201, nlminb iterations: 7
reStruct  parameters:
       ID1        ID2        ID3 
0.96243507 0.07618717 1.63736945 
PNLS step: RSS =  63.22991 
 fixed effects: -3.211823  0.4451494  -0.7862432  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
      fixed    reStruct 
0.001570583 0.013884126 

**Iteration 5
LME step: Loglik: -179.3277, nlminb iterations: 4
reStruct  parameters:
       ID1        ID2        ID3 
0.96197111 0.07396744 1.63796900 
PNLS step: RSS =  63.24112 
 fixed effects: -3.211722  0.4454457  -0.7862145  
 iterations: 6 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.0006651255 0.0147544693 

**Iteration 6
LME step: Loglik: -179.3246, nlminb iterations: 1
reStruct  parameters:
       ID1        ID2        ID3 
0.96215433 0.07487426 1.63772430 
PNLS step: RSS =  63.24112 
 fixed effects: -3.211722  0.4454457  -0.7862145  
 iterations: 1 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.000000e+00 3.910435e-10 
Warning message:
In nlmixrUI.nlme.var(obj) :
  Initial condition for additive error ignored with nlme
> fit.nlme.ode <- nlmixr(uif.ode, theo_sd, est="nlme", calc.resid=FALSE)

**Iteration 1
LME step: Loglik: -182.2318, nlminb iterations: 1
reStruct  parameters:
      ID1       ID2       ID3 
1.0022618 0.2783008 2.0759012 
PNLS step: RSS =  63.24481 
 fixed effects: -3.211675  0.4451503  -0.7862385  
 iterations: 7 
Convergence crit. (must all become <= tolerance = 1e-05):
    fixed  reStruct 
0.2718787 2.7010441 

**Iteration 2
LME step: Loglik: -179.3267, nlminb iterations: 6
reStruct  parameters:
       ID1        ID2        ID3 
0.96197025 0.07396818 1.63806766 
PNLS step: RSS =  63.23773 
 fixed effects: -3.211876  0.4453302  -0.7862234  
 iterations: 2 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.0004040574 0.0088556330 

**Iteration 3
LME step: Loglik: -179.3253, nlminb iterations: 1
reStruct  parameters:
       ID1        ID2        ID3 
0.96208503 0.07453515 1.63782577 
PNLS step: RSS =  63.23773 
 fixed effects: -3.211876  0.4453302  -0.7862234  
 iterations: 1 
Convergence crit. (must all become <= tolerance = 1e-05):
       fixed     reStruct 
0.000000e+00 3.183879e-10 
Warning message:
In nlmixrUI.nlme.var(obj) :
  Initial condition for additive error ignored with nlme
> # SAEM fittings
> 
> # For SAEM, it's nice to avoid printing each C++ iteration to the monitor.
> sink("saem output")
> fit.saem.lin <- nlmixr(uif.lin, theo_sd, est="saem")  # works with spaces in path
Compiling SAEM user function...C:/RTools/3.4/mingw_64/bin/g++  -I"C:/R/R-34~1.3/include" -DNDEBUG       -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/nlmixr/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/STANHE~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/Rcpp/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/RCPPAR~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/RCPPEI~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/BH/include   -O2 -Wall  -mtune=generic -c saem137829e957ffx64.cpp -o saem137829e957ffx64.o
C:/RTools/3.4/mingw_64/bin/g++ -shared -s -static-libgcc -o saem137829e957ffx64.dll tmp.def saem137829e957ffx64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
done.
Using sympy via SnakeCharmR
C:/RTools/3.4/mingw_64/bin/gcc  -I"C:/R/R-34~1.3/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_920cdb7f8c8c7b5e42b1a8db66863398_x64.c -o rx_920cdb7f8c8c7b5e42b1a8db66863398_x64.o
C:/RTools/3.4/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_920cdb7f8c8c7b5e42b1a8db66863398_x64.dll tmp.def rx_920cdb7f8c8c7b5e42b1a8db66863398_x64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
The model has been setup to calculate residuals.
It will be cached for future runs.
Calculating Table Variables...
done
nlmixr UI combined dataset and properties
 $ par.hist         : Parameter history (if available)
 $ par.hist.stacked : Parameter history in stacked form for easy plotting (if available)
 $ par.fixed        : Fixed Effect Parameter Table
 $ eta              : Individual Parameter Estimates
 $ seed             : Seed (if applicable)
 $ model.name       : Model name (from R function)
 $ data.name        : Name of R data input
> # Using RStudio, in Environment tab, what is object curi with value 0 has appeared.
> # What is curi?
> # NB: The following won't work with spaces in file full file path, don't know why:
> fit.saem.ode <- nlmixr(uif.ode, theo_sd, est="saem")  # fails with spaces in path
Compiling RxODE differential equations...done.
C:/RTools/3.4/mingw_64/bin/g++  -I"C:/R/R-34~1.3/include" -DNDEBUG       -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/nlmixr/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/STANHE~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/Rcpp/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/RCPPAR~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/RCPPEI~1/include -IC:/Users/James/DOCUME~1/R/WIN-LI~1/3.4/BH/include   -O2 -Wall  -mtune=generic -c saem13783d7445b6x64.cpp -o saem13783d7445b6x64.o
C:/RTools/3.4/mingw_64/bin/g++ -shared -s -static-libgcc -o saem13783d7445b6x64.dll tmp.def saem13783d7445b6x64.o C:/Users/James/Documents/pharmacometrics/practice/nlmixr/rx_7cbddd1a680b22377f3ee375ef2ddda8_x64.dll -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
done.
## Calculate ETA-based prediction and error derivatives:
Calculate Jacobian...................done.
Calculate sensitivities.......
done.
## Calculate d(f)/d(eta) 
## ...
## done
## ...
## done
C:/RTools/3.4/mingw_64/bin/gcc  -I"C:/R/R-34~1.3/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_32a3413ee26a7ebdb0539062d9fc430e_x64.c -o rx_32a3413ee26a7ebdb0539062d9fc430e_x64.o
C:/RTools/3.4/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_32a3413ee26a7ebdb0539062d9fc430e_x64.dll tmp.def rx_32a3413ee26a7ebdb0539062d9fc430e_x64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
C:/RTools/3.4/mingw_64/bin/gcc  -I"C:/R/R-34~1.3/include" -DNDEBUG          -O2 -Wall  -std=gnu99 -mtune=generic -c rx_a2baf9a7f76ed6696b7c9c5507ef49d7_x64.c -o rx_a2baf9a7f76ed6696b7c9c5507ef49d7_x64.o
C:/RTools/3.4/mingw_64/bin/gcc -shared -s -static-libgcc -o rx_a2baf9a7f76ed6696b7c9c5507ef49d7_x64.dll tmp.def rx_a2baf9a7f76ed6696b7c9c5507ef49d7_x64.o -LC:/R/R-34~1.3/bin/x64 -lRblas -LC:/R/R-34~1.3/bin/x64 -lRlapack -lgfortran -lm -lquadmath -LC:/R/R-34~1.3/bin/x64 -lR
The model-based sensitivities have been calculated.
It will be cached for future runs.
Calculating Table Variables...
done
nlmixr UI combined dataset and properties
 $ par.hist         : Parameter history (if available)
 $ par.hist.stacked : Parameter history in stacked form for easy plotting (if available)
 $ par.fixed        : Fixed Effect Parameter Table
 $ eta              : Individual Parameter Estimates
 $ seed             : Seed (if applicable)
 $ model.name       : Model name (from R function)
 $ data.name        : Name of R data input
> sink()
> # Now compare fits for explicit (i.e., linCmt() ) model vs. ODE, and by algorithms:
> fit.nlme.lin  # Log-likelihood: -179.3246
Nonlinear mixed-effects model fit by maximum likelihood
  Model: DV ~ (nlmeModList("user_fn"))(logtcl, eta.cl, logtka, eta.ka,      logtv, eta.v, TIME, ID) 
  Log-likelihood: -179.3246
  Fixed: logtcl + logtka + logtv ~ 1 
    logtcl     logtka      logtv 
-3.2117217  0.4454457 -0.7862145 

Random effects:
 Formula: list(eta.cl ~ 1, eta.ka ~ 1, eta.v ~ 1)
 Level: ID
 Structure: Diagonal
           eta.cl    eta.ka    eta.v  Residual
StdDev: 0.2644566 0.6422369 0.134573 0.6921699

Number of Observations: 132
Number of Groups: 12 
> fit.nlme.ode  # Log-likelihood: -179.3253
Nonlinear mixed-effects model fit by maximum likelihood
  Model: DV ~ (nlmeModList("user_fn"))(logtcl, eta.cl, logtka, eta.ka,      logtv, eta.v, TIME, ID) 
  Log-likelihood: -179.3253
  Fixed: logtcl + logtka + logtv ~ 1 
    logtcl     logtka      logtv 
-3.2118758  0.4453302 -0.7862234 

Random effects:
 Formula: list(eta.cl ~ 1, eta.ka ~ 1, eta.v ~ 1)
 Level: ID
 Structure: Diagonal
           eta.cl    eta.ka     eta.v  Residual
StdDev: 0.2644679 0.6424376 0.1345558 0.6921515

Number of Observations: 132
Number of Groups: 12 
> fit.saem.lin  # Log-likelihood: -148.4444
nlmixr SAEM fit (Solved)

     OBJF      AIC      BIC Log-likelihood
 296.8482 310.8482 331.0278      -148.4241

Time (sec; $time):
        saem setup Likelihood Calculation covariance table
elapsed 53.9 65.11                   2.09          0   1.5

Parameters ($par.fixed):
                          Parameter Estimate     SE    CV Untransformed
logtcl  log typical value Cl (L/hr)    -3.22 0.0816  2.5%        0.0401
logtka  log typical value Ka (1/hr)    0.450  0.194 43.0%          1.57
logtv      log typical value V  (L)   -0.784 0.0435  5.6%         0.457
add.err                     add.err    0.692                      0.692
                 (95%CI)
logtcl  (0.0342, 0.0471)
logtka      (1.07, 2.29)
logtv     (0.419, 0.497)
add.err                 

Omega ($omega):
          eta.cl    eta.ka      eta.v
eta.cl 0.0696184 0.0000000 0.00000000
eta.ka 0.0000000 0.4187747 0.00000000
eta.v  0.0000000 0.0000000 0.01819676

Fit Data (object is a modified data.frame):
# A tibble: 132 x 15
  ID     TIME    DV IPRED  PRED   IRES    RES  IWRES     WRES  CWRES   CPRED     CRES
  <fct> <dbl> <dbl> <dbl> <dbl>  <dbl>  <dbl>  <dbl>    <dbl>  <dbl>   <dbl>    <dbl>
1 1     0     0.740  0     0     0.740 0.740   1.07   2.06e-5 -0.999  4.81e4 -4.81e+4
2 1     0.250 2.84   3.85  2.82 -1.01  0.0178 -1.46   2.57e-2 -1.46   3.85e0 -1.01e+0
3 1     0.570 6.57   6.76  5.06 -0.191 1.51   -0.276  2.19e+0 -0.276  6.76e0 -1.91e-1
# ... with 129 more rows, and 3 more variables: eta.cl <dbl>, eta.ka <dbl>,
#   eta.v <dbl>
> fit.saem.ode  # Log-likelihood:  -58.07448
nlmixr SAEM fit (ODE)

    OBJF     AIC      BIC Log-likelihood
 116.149 130.149 150.3286      -58.07448

Time (sec; $time):
         saem setup Likelihood Calculation covariance table
elapsed 49.99 17.71                   0.28          0  1.18

Parameters ($par.fixed):
                          Parameter Estimate     SE    CV Untransformed
logtcl  log typical value Cl (L/hr)    -3.22 0.0818  2.5%        0.0400
logtka  log typical value Ka (1/hr)    0.458  0.192 41.9%          1.58
logtv      log typical value V  (L)   -0.782 0.0437  5.6%         0.458
add.err                     add.err    0.694                      0.694
                 (95%CI)
logtcl  (0.0341, 0.0470)
logtka      (1.09, 2.30)
logtv     (0.420, 0.499)
add.err                 

Omega ($omega):
           eta.cl    eta.ka      eta.v
eta.cl 0.07190448 0.0000000 0.00000000
eta.ka 0.00000000 0.4174639 0.00000000
eta.v  0.00000000 0.0000000 0.01818897

Fit Data (object is a modified data.frame):
# A tibble: 132 x 15
  ID     TIME    DV IPRED  PRED   IRES     RES  IWRES    WRES  CWRES CPRED  CRES
  <fct> <dbl> <dbl> <dbl> <dbl>  <dbl>   <dbl>  <dbl>   <dbl>  <dbl> <dbl> <dbl>
1 1     0     0.740  0     0     0.740 0.740    1.07  1.07    1.07    0    0.740
2 1     0.250 2.84   3.90  2.83 -1.06  0.00644 -1.53  0.00382 0.0813  2.66 0.177
3 1     0.570 6.57   6.83  5.07 -0.255 1.50    -0.368 0.675   0.632   4.83 1.74 
# ... with 129 more rows, and 3 more variables: eta.cl <dbl>, eta.ka <dbl>,
#   eta.v <dbl>
> # Why are they so very different?
> # Why does object x appear (1 obs. of 4 variables) appear after fit.saem.ode only?
> 
> # Suggestion: Use same printing format for models with either algorithm.`

nlme method leaves temporary folders behind

The nlme_ode method leaves its temporary folders behind (e.g. file1b9d450716952.d) in the working folder. These should ideally be tidied up after estimation has been concluded. Also files fn.txt, Makevars. If we want to keep these temp folders, can we provide user option to name them and/or autoname them something non-random?

DV in multiple compartments

Dear nlmixr-team

In a permeability assay I have observed concentrations in two compartments (either side of the mebrane).
I included a simple structural model with error structure for both compartments as indicated below.

model({
    Q12 <- exp(tQ12)
    Q21 <- exp(tQ21)
    d/dt(Cap) = -Q12 * Cap + Q21 * Cbl
    d/dt(Cbl) =  Q12 * Cap - Q21 * Cbl
    Cap ~ prop(prop.err)
    Cbl ~ prop(prop.err)
  })

But following errors occur:

  • with "nlme" : Error in varPower(form = ~fitted(.) | nlmixr.grp, fixed = c(1)) : fixed parameters must have group names in 'varPower'
  • with "saem" : Error in nlmixrUI.saem.res.mod(obj) : Currently SAEM only supports one type of residual error model.

Is it currently possible to use nlmixr with observations in 2 compartments? Or is this related to something else? Any additional thoughts?
Thx!

Unclear how to specify covariates

It is unclear how to specify covariates in the ODE form.

The following code works perfectly fine:

# NLMixr model-based analysis; no covariates
library(nlmixr, quietly=TRUE)
dat <- read.table(system.file("examples/theo_md.txt", package = "nlmixr"), head=TRUE)
ode <- "
KE = CL/V;
d/dt(depot) =-KA*depot;
d/dt(centr) = KA*depot - KE*centr;
CONC = centr/V;
"
mypar <- function(lKA, lV, lCL)
{
  KA <- exp(lKA)
  V <- exp(lV)
  CL <- exp(lCL)
}
specs <- list(fixed=lKA+lV+lCL~1, 
              random = pdDiag(lKA+lCL~1), 
              start=c(lKA=0.5, lV=-2.5, lCL=-3.2))
fit <- nlme_ode(
  dat, 
  model=ode, 
  par_model=specs, 
  par_trans=mypar, 
  response="CONC"
)
nlme_gof(fit)

Adding in a covariate:

# NLMixr model-based analysis; covariates
library(nlmixr, quietly=TRUE)
dat <- read.table(system.file("examples/theo_md.txt", package = "nlmixr"), head=TRUE)
ode <- "
KE = (WT/70)^0.75 * CL/V;
d/dt(depot) =-KA*depot;
d/dt(centr) = KA*depot - KE*centr;
CONC = centr/V;
"
mypar <- function(lKA, lV, lCL)
{
  KA <- exp(lKA)
  V <- exp(lV)
  CL <- exp(lCL)
}
specs <- list(fixed=lKA+lV+lCL~1, 
              random = pdDiag(lKA+lCL~1), 
              start=c(lKA=0.5, lV=-2.5, lCL=-3.2))
fit <- nlme_ode(
  dat, 
  model=ode, 
  par_model=specs, 
  par_trans=mypar, 
  response="CONC"
)
nlme_gof(fit)

We then receive an error message:

Error in ..ModList$m1$run(theta, ev, inits, transit_abs = FALSE, atol = 1e-08,  : 
  var(s) not found in input pars.
 WT 

Uniform interface for estimation methods

It would be a huge improvement if it were possible to provide a unified user interface for the different estimation algorithms. Naturally, as far as possible, because the different algorithms may build on different R packages: nlme for instance has its very well developed syntax that you would want to keep being accessible, but some harmonisation should definitely be possible.

For instance: nlme requires parameters to be estimated in the log domain, a natural choice because it allows IIV to be implemented as additive on the log-scale corresponding to exponential-type IIV estimates. And while the SAEM module does basically the same, here parameter are specified on the original scale and are log-transformed internally.

nlme closed-form:

specs1<-list(fixed=lCL+lV~1, random = pdDiag(lCL+lV~1), start=c(lCL=1.6,lV=4.5))
fit <- nlme_lin_cmpt(dat, par_model=specs1, ncmt=1, verbose=TRUE,oral=FALSE,weight=varPower(fixed=c(1)))

nlme-ODE:

ode1 <- "
d/dt(centr)  = -(CL/V)*centr;
"
mypar1 <- function(lCL, lV )
{
    CL <- exp(lCL) 
    V <- exp(lV)
}
fitODE <- nlme_ode(dat, model=ode1, par_model=specs1, par_trans=mypar1, response="centr", response.scaler="V", verbose=TRUE,weight=varPower(fixed=c(1)))

SAEM closed-form:

saem_fit <- gen_saem_user_fn(model=lincmt(ncmt=1, oral=FALSE))
model = list(saem_mod=saem_fit, res.mod=2,covars="WT")
inits = list(theta=c(5,90),omega=c(0.1,0.1),bres=0.2)
cfg   = configsaem(model, nmdat, inits)
cfg$print = 50
fit = saem_fit(cfg)

SAEM ODE:

ode1 <- "
d/dt(centr)  = -(CL/V)*centr;
"
m1 = RxODE(ode1, modName="m1")
mypars1 <- function(lCL, lV )
{
    CL = exp(lCL) 
    V = exp(lV)
}
PRED = function() centr / V
dyn.unload("saem_main.dll")
saem_fit <- gen_saem_user_fn(model=m1, PKpars=mypars1, pred=PRED)

A suggestion would be:
-syntax to generate a compiled model with associated handle (e.g. m1)
-standard syntax to define parameter transformations
-standard syntax to provide inter-individual variability models
-standard syntax to describe residual variability models
-standard syntax to provide starting values
-standard syntax to start a model fit

Would it be an option to use exisiting nlme model definition as a starting point for this? It is well established syntax and it is probably easier to adapt new implementations to the exisiting standard than to try to adapt nlme syntax to newly developed conventions...

Issues with installation instructions under OSX

In the Installing_nlmixr.pdf file they are a couple issues for the OSX instructions:

  1. Is it currently not possible to copy and paste commands in R from the online pdf version on github > Suggest to change to an HTML format
  2. To install packages from github I believe that the syntax should be: devtools::install_github('nlmixrdevelopment/nlmixr ') not devtools::install_github('nlmixrdevelopment/nlmixr', repos = NULL, type = 'source') which gave erros. The same thing applies to the RxODE.
  3. For RxODE I got the following error even though I have gfortran (v.6.1) installed.
* installing *source* package ‘RxODE’ ...
** libs
** arch - 
echo "make libodeaux.a in ode/ ..."
make libodeaux.a in ode/ ...
(cd ode; make CC="`"/Library/Frameworks/R.framework/Resources/bin/R" CMD config CC`" FC="`"/Library/Frameworks/R.framework/Resources/bin/R" CMD config FC`" CFLAGS="`"/Library/Frameworks/R.framework/Resources/bin/R" CMD config CFLAGS` -fPIC" FFLAGS="`"/Library/Frameworks/R.framework/Resources/bin/R" CMD config FFLAGS` -fPIC" libodeaux.a)
gfortran-4.8 -g -O2 -fPIC  -c -o opkda2.o opkda2.f
make[1]: gfortran-4.8: No such file or directory
make[1]: *** [opkda2.o] Error 1
make: *** [libodeaux] Error 2
ERROR: compilation failed for package ‘RxODE’

I was able to solve the issue using this StackOverflow post.

  1. I just came to realize that RxODE is now available on CRAN and this should be mentioned in the instructions it would prevent the installation issues.

Print similar fit information when residuals are not calculated (was print.saemFit suggestions)

I think the summary output from saemFit can be optimized. Example:

THETA:
             th    log(th) se(log_th)
[1,] 0.01043238 -4.5628413 0.02885907
[2,] 3.50889395  1.2553009 0.02954748
[3,] 0.74304072 -0.2970044 0.03018968
[4,] 0.61105778 -0.4925638 0.11929261
[5,] 0.11811888 -2.1360637 0.40655670

OMEGA:
            [,1]       [,2]       [,3]      [,4]    [,5]
[1,] 0.007456519 0.00000000 0.00000000 0.0000000 0.00000
[2,] 0.000000000 0.02114364 0.00000000 0.0000000 0.00000
[3,] 0.000000000 0.00000000 0.01872651 0.0000000 0.00000
[4,] 0.000000000 0.00000000 0.00000000 0.1291693 0.00000
[5,] 0.000000000 0.00000000 0.00000000 0.0000000 1.55202

ARES & BRES:
          [,1]
[1,] 0.3172321
[2,] 0.1543209
  • it is not obvious to me what ARES and BRES mean (is it proportional and residual?)
  • currently I don't know which parameters are linked to which theta's / omega's. Is it the order of mention in the model? Would prefer to have a column with the parameter name.
  • would be nice to have a table which shows the IIV per parameter, e.g. add the %IIV (and perhaps 95% parametric CI) to the THETA table.

mrgsolve and nlmixr

I have been attempting to write the ODE's and model inputs using mrgsolve while using nlmixr for parameter estimation, similar to this link: metrumresearchgroup/mrgsolve#307

However, Im working with a two-compartment ODE-based model and prefer to use saem over nlme. Saem() requires the ODE's be in RxODE format, so I think it would be easier to use nlmixr(est="saem") as it allows additional model classes to be used. I was wondering if you had any examples similar to the above link, but using the saem method instead.

Option to suppress console compilation output

For automated qualification/testing purposes it would be useful to be able to suppress the R console output of the various fitting functions, in order to generate a clean test log. Can we have flags to enable this?

Functional tests for nlmixr and question about one example, saem_theo1.R

I'm working to install nlmixr for a client of mine, and I'm not sure if the installation was successful. Is there a set of functional tests for nlmixr that I could run to confirm the installation is correct?

Also, I've tried running one of the examples, saem_theo1.R, and it appears to be broken. I get this error:

Error in numeric(neq) : invalid 'length' argument

in the statement, cfg = configsaem(model, nmdat, inits)

This statement is problematic:

neq = attr(model$saem_mod, "neq")

I'm using the version on CRAN running on R-3.4.3 on Red Hat Linux 6.9.

Thank you.

access to internal variables?

was wondering if there is a way in the model definition to use internal variables, especially time (either discrete or continuous time, i.e. TIME or T in NONMEM), and the event type at a specific time? Use case is e.g. calculation of time after dose and interpolation of time-varying covariates. If there are other ways to implement those I'd be happy to learn as well :)

Was also wondering how nlmixr/rxode identifies what part of model definition goes into $PK and what goes into $DES (to use the NONMEM analogy), or does it treat all the code as $DES-type?

Thanks!

Object x sometimes appears after fitting, sometimes not

Hi Matt,

At the end of Issue #29 you wrote that you couldn't get the small object x to appear. (x is a 1 x 1 matrix.) Instead of continuing that issue, which you marked as closed, I thought it would be better to open a new one for it. To get the object to appear, start with a new folder that is empty, do a fit, close RStudio session, and then go back to that same folder (which now won't be empty) and repeat this. The first time you should see object x appear but the second time you shouldn't. I attach my results that bear this out.

James
first_output.txt
second_output.txt

Cannot run two SAEM jobs simultaneously

In contrast to nlmixr/nlme it is currently not possible to run two or more nlmixr/SAEM instances at the same time under windows. For instance, when using Microsoft R Open with nlmixr, there is the possibility of setting up a temporary cluster and running nlme jobs in parallel, like with the following syntax:

library(doParallel)
library(foreach)

cl <- makeCluster(15)
registerDoParallel(cl)

do_nlmixr<-function(i){
datr <- read.csv(paste("D:\\nmrun\\nmrun4\\runN024\\m1\\bs_pr1_",i,".dta",sep=""), header=TRUE,stringsAsFactors=F)
fit <- nlme_lin_cmpt(datr, par_model=specs4b, ncmt=1, verbose=TRUE,oral=TRUE,weight=varPower(fixed=c(1)),control = nlmeControl(pnlsTol = .1))
save(fit,file=paste("fit_",i,".Rdata",sep=""))
}

nlmixr_out<-foreach(i=1:500, .combine='rbind', .packages=c('RxODE','nlmixr','nlme')) %dopar% do_nlmixr(i)

When you try to do something like that with SAEM, it crashes immediately, because a single dll is generated in the run directory that is then shared (?): nlme runs create temporary directories that store the material, but SAEM does not. It would be extremely useful if these processes were separated.

Fitting sparse PK data

Dear All,

nlmixr works pretty well with rich PK data. I wonder what is performance with relative sparse sampling schemes. Like using codes below:

library(nlmixr)

f <- function(){
ini({
lCl <- -3 # log Cl (L/hr)
lVc <- -1 # log Vc (L)
lKA <- 1 # log Ka (1/hr)
add.err <- c(0, 0.2, 1)
eta.Cl ~ 0.1 # BSV Cl
eta.Vc ~ 0.1 # BSV Vc
eta.KA ~ 0.1 # BSV Ka
})
model({
## First parameters are defined in terms of the initial estimates
## parameter names.
Cl <- exp(lCl + eta.Cl)
Vc <- exp(lVc + eta.Vc)
KA <- exp(lKA + eta.KA)
## Next, the differential equations are defined
d/dt(depot) = -KAdepot;
d/dt(centr) = KA
depot-Cl/Vc*centr;
## And the concentration is then calculated
cp = centr / Vc;
## Last, nlmixr is told that the plasma concentration follows
## an additive error (estimated by the parameter add.err)
cp ~ add(add.err)
})
}

fit <- nlmixr(f, theo_sd, est="saem")
plot(fit)
A<-data.frame(theo_sd)
B<-subset(A,TIME>1)
B<-subset(B,TIME<4|TIME>12)
fit1 <- nlmixr(f, B, est="saem")

It doesn't fit for me, is there any other function or coding I can try?

Thanks,
Wangda

Nonlinear non-population ODE model

Hello nlmixr development team,

I am trying out nlmixr for the first time and have an issue. I have a multi dose data set (simulated) for a single subject. I am fixing most of the parameters and trying to estimate one parameter only. I do not need to BSV random effect. I am using nlme_ode() function and am getting the following error message. I believe it is looking for random effect parameters formula. Any suggestions?

Also, how do I mention residual error structure in this type of model?

I am not sure how to attach a dataset here. I can send it in email if needed.

Thanks,
Krina

Model Code:

ode <- "
CL=0.111;
Vmax=11.9;
Km=33.9;
V1=2.91;
Q=0.445;
V2=3.06;
C1 = centr/V1;
C2 = peri/V2;
d/dt(depot) = -KA*depot;
d/dt(centr) = KA*depot - CL*C1 - Vmax*C1/(Km+C1) - Q*C1 + Q*C2;
d/dt(peri)  = Q*C1 - Q*C2;
"
specs <- list(fixed=lKA~1, start=c(KA=-0.916))

mypar <- function(lKA) {KA=exp(lKA)}

fit <- nlme_ode(fitdata, ode, par_model=specs, par_trans=mypar, response="C1")

Error in reStruct(random, REML = REML, data = NULL) :
'object' must be a list or a formula

If I update the specs with the following, it works
specs <- list(fixed=lKA~1, random=lKA~1, start=c(KA=-0.916))

Error running nlmixr

Hello,

I've been trying to run the example to test the newly installed nlmixr package but keeping seeing this error message quickly after this part:

fit <- nlmixr(uif, theo_sd, est="nlme", calc.resid=FALSE)
Error: 'rxIs' is not an exported object from 'namespace:RxODE'

I've even tried reinstalling but still no luck. Does anyone know what's the problem?

miscellaneous issues fitting

Hi all, a few more things I encountered during testing. Instead of separating into different issues I combined them in single document, as I'm not sure if they're somehow related or not.

See docs (md + csv) here: https://gist.github.com/ronkeizer/8d34453fc3b133a570ffc56581519146#file-test-rmd

In summary:

  • SAEM works great, but calc.residuals not working / taking very long time
  • FOCEI seems unable to handle infusions events
  • nlme and focei both crash on the full model. A reduced model does slightly better, but still crashes after many iterations.

Issues with "Running PK models with nlmixr" tutorial

Hi Matt,

As I think my nlmixr is properly running as far as I can tell (as I just posted in Issue #25), I think that there are two other issues with the "Running PK models with nlmixr" tutorial at https://cran.r-project.org/web/packages/nlmixr/vignettes/running_nlmixr.html which need addressing. Following along with the first fit with the nlme algorithm, the numbers are different (between what I get and what's in the tutorial) but the differences are small, so that may be a bit irksome but isn't too disconcerting. However, with the SAEM algorithm, the differences are substantial, as I reported a few weeks ago (but maybe the emphasis was lost with other things). For example, in the tutorial the objective function value was 116.211 and the log likelihood was -58.1055, whereas I got OBJF = 349.6266 and log-likelihood = -174.8133 respectively. I did set.seed() to two different values but seem to have gotten the same results both times. Also, the plot should be unimodal for the theophylline data, but it's bimodal because somehow other data were included at a much later time. (The theophylline data should only go up to about 25 hours, not 150.) Here's the pasted output for the plot; there is a problem here somehow—DV shouldn't be lost:

plot(fit.saem)
geom_smooth() using method = 'loess'
Warning messages:
1: Computation failed in stat_smooth():
object 'DV' not found
2: Computation failed in stat_smooth():
object 'DV' not found

tutorial output.txt
Rplot should be unimodal for fit.saem.pdf

CRAN version does not work with nlme estimation

Hi!

Saw NLMIXR was now on CRAN. Tested the source version on Linux (R3.4.1) and Windows and the binary version on Windows in R3.3.3 and 3.4.1.

None of the versions worked on my machines and on my example. I had a self compiled version previously which worked nicely on all my NLMIXR examples ... bit weird ... especially as it had the same version number as the one on CRAN now.

Example available here: http://web.intiquan.com/Stuff/NLMIXR/exampleNLMIXR.zip

Cheers,
Henning

SymPy Installation Issues

I installed and updated all recommended packages from the installation page on a MacOS Sierra. When attempting to run nlmixr examples, there is an issue loading sympy via snakecharmR . I keep getting the following error message:

Error in SnakeCharmR::py.exec(...) : Traceback (most recent call last):
File "", line 3, in
ImportError: No module named sympy

/usr/bin/python: No module named pip
Error in SnakeCharmR::py.exec(...) : Traceback (most recent call last):
File "", line 3, in
ImportError: No module named sympy

I'm assuming the error has to do with the Path snakecharmR uses when retrieving sympy. I installed pip, SnakeCharmR, and sympy successfully, but I can't figure out why SnakeCharmR isn't able to load sympy. Any suggestions on how to fix this issue/ update the path?

Error installing SnakeCharmR

Hi,

I am trying to install SnakeCharmR in a RedHat Linux environment.

Command used: install_github("nlmixrdevelopment/SnakeCharmR

Error Log:

* installing *source* package ‘SnakeCharmR’ ...
** package ‘SnakeCharmR’ successfully unpacked and MD5 sums checked
Using python-config...
Traceback (most recent call last):
  File "<string>", line 1, in <module>
AttributeError: 'tuple' object has no attribute 'major'
Found Python version ...
configure: creating ./config.status
config.status: creating src/Makevars
--------[begin src/Makevars]--------
# Required vars
PKG_LIBS=-lpthread -ldl -lutil -lm -lpython2.6
PKG_CPPFLAGS=-I/usr/include/python2.6 -I/usr/include/python2.6 -D PYTHONLIBFILE=\"libpython.so\"
CXX_STD = CXX11
--------[end src/Makevars]--------
** libs
g++ -m64 -std=gnu++0x -I/usr/include/R -DNDEBUG -I/usr/include/python2.6 -I/usr/include/python2.6 -D PYTHONLIBFILE=\"libpython.so\" -I"/      usr/lib64/R/library/Rcpp/include" -I/usr/local/include   -fpic  -O2 -g -pipe -Wall -Wp,-D_FORTIFY_SOURCE=2 -fexceptions -fstack-protecto      r --param=ssp-buffer-size=4 -m64 -mtune=generic -c RcppExports.cpp -o RcppExports.o
g++ -m64 -std=gnu++0x -I/usr/include/R -DNDEBUG -I/usr/include/python2.6 -I/usr/include/python2.6 -D PYTHONLIBFILE=\"libpython.so\" -I"/      usr/lib64/R/library/Rcpp/include" -I/usr/local/include   -fpic  -O2 -g -pipe -Wall -Wp,-D_FORTIFY_SOURCE=2 -fexceptions -fstack-protecto      r --param=ssp-buffer-size=4 -m64 -mtune=generic -c python.cpp -o python.o
In file included from /usr/include/python2.6/pyconfig.h:6:0,
                 from /usr/include/python2.6/Python.h:8,
                 from python.cpp:4:
/usr/include/python2.6/pyconfig-64.h:1034:0: warning: "_POSIX_C_SOURCE" redefined
 #define _POSIX_C_SOURCE 200112L
 ^
In file included from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/x86_64-redhat-linux/bits/os_defines.h:39:0,
                 from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/x86_64-redhat-linux/bits/c++config.h:2255,
                 from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/cmath:41,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp/platform/compiler.h:100,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp/r/headers.h:48,
                 from /usr/lib64/R/library/Rcpp/include/RcppCommon.h:38,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp.h:27,
                 from python.cpp:1:
/usr/include/features.h:162:0: note: this is the location of the previous definition
 # define _POSIX_C_SOURCE 200809L
 ^
In file included from /usr/include/python2.6/pyconfig.h:6:0,
                 from /usr/include/python2.6/Python.h:8,
                 from python.cpp:4:
/usr/include/python2.6/pyconfig-64.h:1043:0: warning: "_XOPEN_SOURCE" redefined
 #define _XOPEN_SOURCE 600
 ^
In file included from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/x86_64-redhat-linux/bits/os_defines.h:39:0,
                 from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/x86_64-redhat-linux/bits/c++config.h:2255,
                 from /opt/rh/devtoolset-4/root/usr/include/c++/5.3.1/cmath:41,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp/platform/compiler.h:100,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp/r/headers.h:48,
                 from /usr/lib64/R/library/Rcpp/include/RcppCommon.h:38,
                 from /usr/lib64/R/library/Rcpp/include/Rcpp.h:27,
                 from python.cpp:1:
/usr/include/features.h:164:0: note: this is the location of the previous definition
 # define _XOPEN_SOURCE 700
 ^
g++ -m64 -std=gnu++0x -shared -L/usr/lib64/R/lib -o SnakeCharmR.so RcppExports.o python.o -lpthread -ldl -lutil -lm -lpython2.6 -L/usr/l      ib64/R/lib -lR
installing to /usr/lib64/R/library/SnakeCharmR/libs
** R
** inst
** preparing package for lazy loading
** help
*** installing help indices
  converting help for package ‘SnakeCharmR’
    finding HTML links ... done
    SnakeCharmR                             html
    py.assign                               html
    py.call                                 html
    py.exec                                 html
    py.get                                  html
    py.load                                 html
    py.rm                                   html
** building package indices
** testing if installed package can be loaded
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "/usr/lib64/python2.6/json/__init__.py", line 109, in <module>
    from .encoder import JSONEncoder
  File "/usr/lib64/python2.6/json/encoder.py", line 5, in <module>
    import math
ImportError: /usr/lib64/python2.6/lib-dynload/mathmodule.so: undefined symbol: PyExc_ValueError
Error: package or namespace load failed for ‘SnakeCharmR’:
 .onAttach failed in attachNamespace() for 'SnakeCharmR', details:
  call: py.get("sys.version")
  error: Traceback (most recent call last):
  File "<string>", line 3, in <module>
NameError: name 'json' is not defined

Error: loading failed
Execution halted
ERROR: loading failed
* removing ‘/usr/lib64/R/library/SnakeCharmR’

The downloaded source packages are in
        ‘/tmp/RtmpKKpNmy/downloaded_packages’
Updating HTML index of packages in '.Library'
Making 'packages.html' ... done
Warning message:
In install.packages("SnakeCharmR") :
  installation of package ‘SnakeCharmR’ had non-zero exit status

Non-Population Single Compartment

Hello,
I'm trying to use this package for modeling non-population data for a single-compartment model. I thought dynmodel might be the function that could help me to that end, but I keep getting an optimHess error when running it.

Does anyone have any guidance on how to use this package for non-population data, or using single-compartment models on the average of non-population data?

Best,
Christopher Neil

Recommend Projects

  • React photo React

    A declarative, efficient, and flexible JavaScript library for building user interfaces.

  • Vue.js photo Vue.js

    🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.

  • Typescript photo Typescript

    TypeScript is a superset of JavaScript that compiles to clean JavaScript output.

  • TensorFlow photo TensorFlow

    An Open Source Machine Learning Framework for Everyone

  • Django photo Django

    The Web framework for perfectionists with deadlines.

  • D3 photo D3

    Bring data to life with SVG, Canvas and HTML. 📊📈🎉

Recommend Topics

  • javascript

    JavaScript (JS) is a lightweight interpreted programming language with first-class functions.

  • web

    Some thing interesting about web. New door for the world.

  • server

    A server is a program made to process requests and deliver data to clients.

  • Machine learning

    Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.

  • Game

    Some thing interesting about game, make everyone happy.

Recommend Org

  • Facebook photo Facebook

    We are working to build community through open source technology. NB: members must have two-factor auth.

  • Microsoft photo Microsoft

    Open source projects and samples from Microsoft.

  • Google photo Google

    Google ❤️ Open Source for everyone.

  • D3 photo D3

    Data-Driven Documents codes.