nonmem summary fails to produce standard error output when not all estimation methods have a $COV step

[bmreilly]

I have a model using SAEM followed by IMP estimation. The IMP estimation had option NOCOV=1 which causes it to skip the covariance step for IMP, but it still produces a covariance matrix for SAEM. There are good reasons to do this: (1) SAEM already computed the covariance matrix so why repeat it; and (2) the $COV step for IMP can take excessively long compared to other methods. This is a fairly common option in NONMEM.

When I call bbr::model_summary() in R or bbi nonmem summary on the command line, neither produces standard errors for the parameter estimates, even though they are available from the SAEM covariance matrix.

Is there some way to tell bbi/bbr to use a specific covariance matrix if there are more than one?

[kyleam]

When I call bbr::model_summary() in R or bbi nonmem summary on the command line, neither produces standard errors for the parameter estimates, even though they are available from the SAEM covariance matrix.

Is there some way to tell bbi/bbr to use a specific covariance matrix if there are more than one?

Both hard code the final estimation method when displaying the result:

• https://github.com/metrumresearchgroup/bbi/blob/835953efed10fff2d7024ac1c3dd88949b268c05/parsers/nmparser/summary.go#L23-L26
• https://github.com/metrumresearchgroup/bbr/blob/1967a101de200e86644f4a3096a408cee9ccf97f/R/param-estimates.R#L64-L72

The values should be there for you to grab with a bit of digging, though, which may be good enough depending on your use case.

with bbi

Default output lacks standard errors, as you describe:

$ bbi nonmem summary 3/3
PK model 1 cmt base
Dataset: ../acop.csv
Records: 779   Observations: 741  Subjects: 39
Estimation Method(s):
 - Stochastic Approximation Expectation-Maximization
 - Objective Function Evaluation by Importance Sampling
No Heuristic Problems Detected
+-------+--------+-----------+--------------+
| THETA |  NAME  | ESTIMATE  | STDERR (RSE) |
+-------+--------+-----------+--------------+
| TH 1  | THETA1 | 2.24689   | -            |
| TH 2  | THETA2 | 47.2305   | -            |
| TH 3  | THETA3 | 15.4904   | -            |
| TH 4  | THETA4 | 0.0549285 | -            |
| TH 5  | THETA5 | 1.34576   | -            |
+-------+--------+-----------+--------------+
+------------+------+-----------+---------------+
|   OMEGA    | ETA  | ESTIMATE  | SHRINKAGE (%) |
+------------+------+-----------+---------------+
| OMEGA(1,1) | ETA1 | 0.136931  | 14.801200     |
| OMEGA(2,2) | ETA2 | 11.741500 | 79.849300     |
+------------+------+-----------+---------------+

But you can extract them from the json output:

$ bbi nonmem summary --json 3/3 | jq '.run_details.estimation_method'
[
  "Stochastic Approximation Expectation-Maximization",
  "Objective Function Evaluation by Importance Sampling"
]

$ bbi nonmem summary --json 3/3 | jq '.parameters_data | .[].std_err'
{
  "theta": [
    0.0664072,
    3.63853,
    3.37799,
    0.00429638,
    0.000541194
  ],
  "omega": [
    0.06093,
    0,
    17.0897
  ],
  "sigma": [
    0
  ]
}
{
  "theta": [
    -999999999,
    -999999999,
    -999999999,
    -999999999,
    -999999999
  ],
  "omega": [
    -999999999,
    -999999999,
    -999999999
  ],
  "sigma": [
    -999999999
  ]
}

with bbr

print.bbi_nonmem_summary doesn't show standard errors:

> bbr::model_summary(mod)
 Dataset: ../acop.csv

 Records: 779     Observations: 741       Subjects: 39

 Objective Function Value (final est. method): 4664.232

 Estimation Method(s):

 – Stochastic Approximation Expectation-Maximization

 – Objective Function Evaluation by Importance Sampling

 **Heuristic Problem(s) Detected:**

 – eta_pval_significant

 # A tibble: 7 × 4
   parameter_names estimate stderr shrinkage
   <chr>           <chr>    <chr>  <chr>
 1 THETA1          2.25     "  "   "  "
 2 THETA2          47.2     "  "   "  "
 3 THETA3          15.5     "  "   "  "
 4 THETA4          0.0549   "  "   "  "
 5 THETA5          1.35     "  "   "  "
 6 OMEGA(1,1)      0.137    "  "   "14.8"
 7 OMEGA(2,2)      11.7     "  "   "79.8"

But you can get them from the summary object:

> s <- bbr::model_summary(mod)
> s[["parameters_data"]][[1]][["std_err"]]
 $theta
 [1] 0.066407200 3.638530000 3.377990000 0.004296380 0.000541194

 $omega
 [1]  0.06093  0.00000 17.08970

 $sigma
 [1] 0

[bmreilly]

@kyleam Thanks for this, this is very helpful. I can see the values in the model summary object as you suggested.

To close the loop on this, I've since learned that the covariance step from SAEM is not as accurate as for other methods, so this workflow should not be very commonly used (although it is expressly mentioned in the NONMEM guide). See: https://www.mail-archive.com/nmusers@globomaxnm.com/msg05861.html

Thanks for your quick and helpful response!