Question regarding extracting parameter uncertainties

[paulrbuckley]

Hi. I am an avid user of dMod, I really like the framework.

I have a model with 20 parameters and I have used dMod and in particular the ms-trust fitting algorithm to derive ~60 sets of parameters which appear to fit the data well. I am trying to generate confidence intervals around each parameter.

I have done this purely from the distribution derived from the multi start fitting algorithm, however I would like to use the profile-likelihood method and the confint function for obvious reasons. The profile-likelihood has shown that many of the model's parameters are not identifiable, however the visualisations through 'plotProfile' show that some cross the threshold in both directions.

For one of these parameters that cross the 95% threshold in both directions, I am still only able to compute a lower confidence bound, with upper returning 'NA'. I am confused about this because visually, it appears to cross the threshold.. This is the case when I 1) compute PL for the best fit parameter set and 2) for the entire parframe.

From delving into some of the code in the function confint.parframe , my data returns NA for the right bound after the following code;

if (i.below > i.above) { return(NA) }

Apologies if i've missed something obvious, but my question is, what does this represent? Further, in such a situation I am wondering whether it is indeed possible to compute the confidence intervals?

Appreciate your comments

cheers

Paul

[paulrbuckley]

Just to add to this, i've just updated dMod from github and I now receive the error "Error in { : task 1 failed - "could not find function "wrss""" when I run the profile function

[marcusrosenblatt]

Dear Paul,

just a quick answer for now: The profile plotting function automatically plots different contributions to your objective function (e.g. data and prior). Make sure that you plot the right contribution by using for example sth like plotProfile(profiles, mode=="data") or plotProfile(profiles, mode=="chisquare")

The wrss function is no longer supported in the current version. Nevertheless the profile function should work properly. Maybe delete your workspace and create all objects again with the new version. Do you then still have the same error?

Best,
Marcus

[paulrbuckley]

Hi Marcus,
Thanks for getting back to me so quickly!

The profile function is working now thank you. I re-ran the entire analysis and that has resolved it.

Regarding the plots for contribution to the objective function, I have tried this just now for both the data and chisq contributions. The threshold is reached in both directions yet I am unable to compute an upper confidence bound for this parameter. I'm a little confused by this!

cheers

P

[marcusrosenblatt]

Hm, could you copy and paste the corresponding parframe. Or even better: attach it as an .RData file or sth?

[paulrbuckley]

Hi,

Sorry for delayed reply! I have managed to compute CIs for the best fit set of parameters, for the parameter mentioned above. This wasn't working on friday but for whatever reason it is working now.

Lastly, I am trying to do this for a parframe (is this possible?) for which I have computed the profiles by:

`profiles <- lapply(1:nrow(parFrameDistribution), function(i) {

pars= as.parvec(parFrameDistribution, i)

profile(obj, pars, names(PAR_WHICH),cores = 6,fixed=c(fixed,pars[!names(pars) %in% names(PAR_WHICH)]),method = "integrate")
})`

where 'names(PAR_WHICH)' evaluates to the name of the parameter in question and 'fixed' are the fixed parameters and initial conditions. I can visualise the profile(s) in this manner, but I get the following error after using 'confint' ;

Error in UseMethod("vcov") : no applicable method for 'vcov' applied to an object of class "list"

Should I convert the 'lapply' resulting list to a data table and then pass it to 'confint', or is there a better way to do this?

cheers