Open Cole-Monnahan-NOAA opened 2 years ago
Another way is to just create another file. Sounds like a good idea. Let me know if you want to proceed.
Yeah we could do it as a separate file I guess, at least for a while and see. It's possible people parse the .std file in externally and those two extra columns would break that code.
I realize now I should have probably put this in as a "discussion" to get more of a ... discussion?
Cole, this is a great plan, thanks for thinking of it. I'm fine with either separate columns or separate file, but I agree that a separate file initially would avoid any unexpected issues with how folks use the std file.
In addition to a new file, do you think it would also be possible to make the new 95% CI values available via a function run within an ADMB model, similar to however the get_covariance_matrix() works in https://github.com/admb-project/admb/blob/main/src/nh99/mod_sd.cpp?
@iantaylor-NOAA could you give more details? In at least the FINAL section you can retrieve the SE of a parameter or derived quanttiy and use that.
dvar_vector UCI = exp(biomsd+1.96*biomsd.sd);
where biomsd is log-biomass in this case, declared as an sdreport vector
Cole, that's bounded SE, right? I'm thinking about how Stock Synthesis reports the bounded SE for all parameters and derived quantities in a custom Report.sso file and it would be nice to also report in that same place the naturally bounded and asymmetric confidence intervals that you propose calculating from the unbounded SE.
admb-13.1 will be released at the before or the beginning of 2023. this would be a nice feature to add in.
Indeed it is in bounded space. All that's needed is the "scale" which is the derivative of the function, in this case the bounding function. So if y=f(x) and x is the unbounded parameter then Var(y)=Var(x)*f'(x)^2. So what we want is to get Var(x) and calculate the CI and then bound that.
The delta method calc is done here.
The header for the .std file (ofsd object) is here.
Then the est and SE values (bounded space) are written on these lines.
So all we need to do, I think, is create a new ofsd object and copy the header, then take the original SE in unbounded space, get CI. Then get the bounds and pass pass the lwr and upr limits to boundp and write them to the new file. Done? Always seems simpler in theory but that's how I'd go about it.
Do development on issue275 branch.
I won't be able to tackle this until December probably, so happy to let someone else give it a shot first.
Okay, I can work on this next month.
I've been thinking recently about how the current approach to calculating confidence intervals is inadequate statistically. The .std file reports the SE in bounded space which can lead to a CI that violates the bounds.
A better approach is to use the unbounded SE and calculate the CI and transform those to the bounded scale. This leads to CIs that are naturally bounded and asymmetric. ADMB has all this information available when the .std file is created (see equations in this document) and so it should be relatively easy to add this to the output.
Specifically, I propose adding two new columns to the .std file: the lower and upper 95% confidence interval limits in bounded space. The user can then use the current approach (backwards compatibility), but also use the new CI in if that is desired. I don't think it's helpful to report the unbounded SE b/c there's no easy way for the user to replicate the bounding function, to, e.g., calculate other coverage intervals, but that of course be reported as well.
Unbounded parameters and derived quantities would be unaffected and the CI would simply be (est-1.96SE, est+1.96SE).
Thoughts?