Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

HRs from a loglog orm do not fall within their corresponding CIs. #56

Closed
yhpua opened this issue Jan 12, 2018 · 4 comments
Closed

HRs from a loglog orm do not fall within their corresponding CIs. #56

yhpua opened this issue Jan 12, 2018 · 4 comments

Comments

@yhpua
Copy link
Contributor

yhpua commented Jan 12, 2018

Dear Prof Harrell,
I've noticed that the hazard ratios from a loglog orm do not fall within their corresponding CIs. Below, I've provided a minimal example.
Thank you for looking into this and for creating the rms package!
Yonghao

create a datadist object

mtcars.dd <- datadist(mtcars)
options(datadist="mtcars.dd")

mod.loglog <- orm (mpg ~ drat + vs , family =loglog, data = mtcars) # loglog CPM

HR ratios are not within their corresponding CIs

HRs for the predictors are probably exp(beta)

summary (mod.loglog)

@yhpua yhpua closed this as completed Jan 12, 2018
@harrelfe harrelfe reopened this Jan 28, 2018
@harrelfe
Copy link
Owner

Sorry to not have gotten to this sooner. The confidence intervals appear to be correct but the hazard ratios are exp(-beta). I'll take a look at the code to try to debug. If you have time to look at it earlier than I get to it that would be great (code for orm and summary.rms).

@yhpua
Copy link
Contributor Author

yhpua commented Feb 2, 2018

Dear Prof. Harrell,

Thank you so much for your reply! I have noticed that the loglog link function in orm.fit.s is given by -log(-log(x )) and not log(-log(x )) as stated in your excellent Regression Modeling Strategies book (Table 15.1, page 362). If the loglog link function is -log(-log(x )), exp(beta) will be given by logP1/logP2 while exp(-beta) will be given by logP2/logP1. If so, does exp(-beta) give the correct Hazard Ratio?

@harrelfe
Copy link
Owner

harrelfe commented Feb 4, 2018

The logic in summary.orm was unnecessarily complex. Thanks for pointing this out. For the PH model (log-log link) I was negating the linear predictor before anti-logging, but was leaving the confidence limits in their original form. So the confidence limits were correct all along, and the point estimate was wrong. I'll have that logic removed for the next submission to CRAN. Changes will be committed to Github momentarily. Let me know if this change does not makes sense.

harrelfe pushed a commit that referenced this issue Feb 4, 2018
…d ratio for orm fit with log-log family (confidence limits never attempted to do this). See #56
@harrelfe harrelfe closed this as completed Feb 6, 2018
@yhpua
Copy link
Contributor Author

yhpua commented Feb 6, 2018

Thank you for making those changes.

clrpackages pushed a commit to clearlinux-pkgs/R-rms that referenced this issue Jan 28, 2019
…dist options and changed default for scat1d.opts

Changes in version 5.1-3 (2019-01-27)
   * lrm: when lrm.fit fails, don't stop altogether but print warning and return fit with fail=TRUE
   * summary.orm: stopped negating linear predictor before computing hazard ratio for orm fit with log-log family (confidence limits never attempted to do this).  See harrelfe/rms#56
   * Design: stop with error message any any design matrix column name is duplicated as with a predictor X1 with factor level 2 (forming name X12 in model.matrix) and a numeric predictor named X12.  See harrelfe/rms#55
   * contrast.rms: for conf.type='simultaneous' shortened fit$coefficients to only have one intercept so will work with multcomp::glht.  Thanks: Matthew Shun-Shin
   * orm: changed label for Unique Y to Distinct Y
   * survplotp: fixed bug in legends
   * predab.resample: handled case where fitter doesn't include fail=FALSE
   * orm.fit: fixed problem with omission of name of fam argument with offset is present, and bug where initial had the wrong length.  Thanks: Tamas Ferenci
   * inst/tests/cph4.r: commented out line that caused failure
   * plotp: moved generic to Hmisc
   * Design: gave better error message re:ordered factors
   * predictrms: remembered to run scored() when producing design matrix
   * Design: allowed for = in value labels in mmnames.  See github.com/harrelfe/rms/issues/29#issuecomment-417901353
   * plot.calibrate.default: added cex.subtitles, riskdist options and changed default for scat1d.opts

(NEWS truncated at 15 lines)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants