search

harrelfe / rms

Regression Modeling Strategies
https://hbiostat.org/R/rms
Other
170 stars 48 forks source link

contrast in probabilities when using lrm #26

Open sdaza opened 8 years ago

sdaza commented 8 years ago

Apparently, I cannot use a function such as 1/(1 + exp(-x)) in the contrast command.

age <- rnorm(200,40,12)
sex <- factor(sample(c('female','male'),200,TRUE))
logit <- (sex=='male') + (age-40)/5
y <- ifelse(runif(200) <= plogis(logit), 1, 0)
f <- lrm(y ~ pol(age,2)*sex)
getProb <- function (x)  { 1/(1 + exp(-x)) }
contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

harrelfe commented 8 years ago

fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Thu, Jun 30, 2016 at 11:48 AM, Sebastian Daza notifications@github.com wrote:

Apparently, I cannot use a function such as 1/(1 + exp(-x)) in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26, or mute the thread https://github.com/notifications/unsubscribe/ABGO2onE3ZVRxd7tx9UUx3URkQebURbuks5qQ_NDgaJpZM4JCW81 .

sdaza commented 8 years ago

Thank you, I got it. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation should be applied, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

harrelfe commented 8 years ago

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza notifications@github.com wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229942984, or mute the thread https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81 .

ruddjm commented 8 years ago

Hi Frank,

I'm puzzled by why using the transformation in print.contrast.rms doesn't do what he wants. Wouldn't doing print(k, fun=plogis) as you suggest return as contrast in probabilities? How is that result different from doing the complex standard error calculations (assuming he did it correctly)?

I copied his original question and your reply below for reference.

Thanks,

JoAnn

Your reply: fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Original question:

Apparently, I cannot use a function such as |1/(1 + exp(-x))| in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 08:53 AM, Frank Harrell wrote:

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza notifications@github.com wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229942984, or mute the thread

https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81 .

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229952667, or mute the thread https://github.com/notifications/unsubscribe/AFJu9vkIeunp1aqUBPEex_EBTYULC5FKks5qRRvWgaJpZM4JCW81.

harrelfe commented 8 years ago

JoAnn,

Why would a nonlinear transformation work in this way? contrast will correctly get you confidence intervals for transformations of contrasts on the linear predictor scale. It has no way to know about differences in nonlinear transformations of the linear predictor. fun= in print.contrast.rms refers to a last-second transformation of the contrast

contrast.rms handles bootstrapped estimates and it may be possible to extend the function to use bootstrap reps for confidence intervals on the right scale.

BUT note that differences in probabilities is not a great scale anyway.

Frank

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 11:39 AM, JoAnn Rudd Alvarez < notifications@github.com> wrote:

Hi Frank,

I'm puzzled by why using the transformation in print.contrast.rms doesn't do what he wants. Wouldn't doing print(k, fun=plogis) as you suggest return as contrast in probabilities? How is that result different from doing the complex standard error calculations (assuming he did it correctly)?

I copied his original question and your reply below for reference.

Thanks,

JoAnn

Your reply: fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Original question:

Apparently, I cannot use a function such as |1/(1 + exp(-x))| in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 08:53 AM, Frank Harrell wrote:

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza <notifications@github.com

wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229942984, or mute the thread

< https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81

.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229952667, or mute the thread < https://github.com/notifications/unsubscribe/AFJu9vkIeunp1aqUBPEex_EBTYULC5FKks5qRRvWgaJpZM4JCW81 .

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229992802, or mute the thread https://github.com/notifications/unsubscribe/ABGO2kGBdYviFLDokCNARXSbdp1LHhHDks5qRULVgaJpZM4JCW81 .

ruddjm commented 8 years ago

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 12:52 PM, Frank Harrell wrote:

JoAnn,

Why would a nonlinear transformation work in this way? contrast will correctly get you confidence intervals for transformations of contrasts on the linear predictor scale. It has no way to know about differences in nonlinear transformations of the linear predictor. fun= in print.contrast.rms refers to a last-second transformation of the contrast Yes, so it will just apply plogis to the betahat and to each of the limits. But what is the difference between this last-second transformation and what Sebastian wants?

contrast.rms handles bootstrapped estimates and it may be possible to extend the function to use bootstrap reps for confidence intervals on the right scale.

BUT note that differences in probabilities is not a great scale anyway. I thought differences on a probability scale are nice because they're more understandable than odds ratios.

By the way, thanks for answering my questions.

Frank

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 11:39 AM, JoAnn Rudd Alvarez < notifications@github.com> wrote:

Hi Frank,

I'm puzzled by why using the transformation in print.contrast.rms doesn't do what he wants. Wouldn't doing print(k, fun=plogis) as you suggest return as contrast in probabilities? How is that result different from doing the complex standard error calculations (assuming he did it correctly)?

I copied his original question and your reply below for reference.

Thanks,

JoAnn

Your reply: fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Original question:

Apparently, I cannot use a function such as |1/(1 + exp(-x))| in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 08:53 AM, Frank Harrell wrote:

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza <notifications@github.com

wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229942984, or mute the thread

<

https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81

.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229952667, or mute the thread <

https://github.com/notifications/unsubscribe/AFJu9vkIeunp1aqUBPEex_EBTYULC5FKks5qRRvWgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229992802, or mute the thread

https://github.com/notifications/unsubscribe/ABGO2kGBdYviFLDokCNARXSbdp1LHhHDks5qRULVgaJpZM4JCW81 .

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230009008, or mute the thread https://github.com/notifications/unsubscribe/AFJu9tLYLeuyG2JRx3JY9H1fH1bPd2woks5qRVPLgaJpZM4JCW81.

sdaza commented 8 years ago

The function is applying plogis to log odds ratios (as @harrelfe pointed out). If the difference is zero, you get a value of 0.5. The difference between predicted probabilities should be also 0.

I agree about the use of odd ratios and logit coefficients. I think we shouldn't use them to communicate our results and models. Look at this.

You can use the package Zelig to simulate predicted values.

harrelfe commented 8 years ago

It applies plogis to the log odds RATIO which doesnt make sense. On Jul 1, 2016 11:00 AM, "JoAnn Rudd Alvarez" notifications@github.com wrote:

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 12:52 PM, Frank Harrell wrote:

JoAnn,

Why would a nonlinear transformation work in this way? contrast will correctly get you confidence intervals for transformations of contrasts on the linear predictor scale. It has no way to know about differences in nonlinear transformations of the linear predictor. fun= in print.contrast.rms refers to a last-second transformation of the contrast Yes, so it will just apply plogis to the betahat and to each of the limits. But what is the difference between this last-second transformation and what Sebastian wants?

contrast.rms handles bootstrapped estimates and it may be possible to extend the function to use bootstrap reps for confidence intervals on the right scale.

BUT note that differences in probabilities is not a great scale anyway. I thought differences on a probability scale are nice because they're more understandable than odds ratios.

By the way, thanks for answering my questions.

Frank

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 11:39 AM, JoAnn Rudd Alvarez < notifications@github.com> wrote:

Hi Frank,

I'm puzzled by why using the transformation in print.contrast.rms doesn't do what he wants. Wouldn't doing print(k, fun=plogis) as you suggest return as contrast in probabilities? How is that result different from doing the complex standard error calculations (assuming he did it correctly)?

I copied his original question and your reply below for reference.

Thanks,

JoAnn

Your reply: fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Original question:

Apparently, I cannot use a function such as |1/(1 + exp(-x))| in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 08:53 AM, Frank Harrell wrote:

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza <notifications@github.com

wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub <https://github.com/harrelfe/rms/issues/26#issuecomment-229942984 , or mute the thread

<

https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81

.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229952667, or mute the thread <

https://github.com/notifications/unsubscribe/AFJu9vkIeunp1aqUBPEex_EBTYULC5FKks5qRRvWgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229992802, or mute the thread

< https://github.com/notifications/unsubscribe/ABGO2kGBdYviFLDokCNARXSbdp1LHhHDks5qRULVgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230009008, or mute the thread < https://github.com/notifications/unsubscribe/AFJu9tLYLeuyG2JRx3JY9H1fH1bPd2woks5qRVPLgaJpZM4JCW81 .

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230010835, or mute the thread https://github.com/notifications/unsubscribe/ABGO2jh56KD4OIMkwvn5bytAGd1RP01kks5qRVXCgaJpZM4JCW81 .

ruddjm commented 8 years ago

Oh, ok, thanks.

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 03:16 PM, Frank Harrell wrote:

It applies plogis to the log odds RATIO which doesnt make sense. On Jul 1, 2016 11:00 AM, "JoAnn Rudd Alvarez" notifications@github.com wrote:

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 12:52 PM, Frank Harrell wrote:

JoAnn,

Why would a nonlinear transformation work in this way? contrast will correctly get you confidence intervals for transformations of contrasts on the linear predictor scale. It has no way to know about differences in nonlinear transformations of the linear predictor. fun= in print.contrast.rms refers to a last-second transformation of the contrast Yes, so it will just apply plogis to the betahat and to each of the limits. But what is the difference between this last-second transformation and what Sebastian wants?

contrast.rms handles bootstrapped estimates and it may be possible to extend the function to use bootstrap reps for confidence intervals on the right scale.

BUT note that differences in probabilities is not a great scale anyway. I thought differences on a probability scale are nice because they're more understandable than odds ratios.

By the way, thanks for answering my questions.

Frank

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 11:39 AM, JoAnn Rudd Alvarez < notifications@github.com> wrote:

Hi Frank,

I'm puzzled by why using the transformation in print.contrast.rms doesn't do what he wants. Wouldn't doing print(k, fun=plogis) as you suggest return as contrast in probabilities? How is that result different from doing the complex standard error calculations (assuming he did it correctly)?

I copied his original question and your reply below for reference.

Thanks,

JoAnn

Your reply: fun is an argument to print.contrast.rms. So do

k <- contrast(...)

print(k, fun=plogis) # or fun=getProb (same thing)

Original question:

Apparently, I cannot use a function such as |1/(1 + exp(-x))| in the contrast command.

age <- rnorm(200,40,12) sex <- factor(sample(c('female','male'),200,TRUE)) logit <- (sex=='male') + (age-40)/5 y <- ifelse(runif(200) <= plogis(logit), 1, 0) f <- lrm(y ~ pol(age,2)*sex) getProb <- function (x) { 1/(1 + exp(-x)) } contrast(f, list(sex='female', age=30), list(sex='male', age=40), fun = getProb)

Any ideas on how to get a contrast in probabilities?

JoAnn Rudd Alvarez http://biostat.mc.vanderbilt.edu/JoAnnAlvarez

On 07/01/2016 08:53 AM, Frank Harrell wrote:

You would have to program complex standard error calculations to be able to do that. Probability differences are highly specific to the adjustment covariate settings so I do not recommend that approach. If you really want to do it, bootstrap the entire process to get confidence intervals.

Frank E Harrell Jr Professor and Chairman School of Medicine

Department of Biostatistics Vanderbilt University

On Fri, Jul 1, 2016 at 8:12 AM, Sebastian Daza <notifications@github.com

wrote:

Thank you, I got. But this is not what I want. I would like to contrast predicted probabilities: if a logit coefficient is negative, the contrast in predicted probabilities should be also negative, this is not what I get. So, the transformation needs to be done, I guess, at the prediction procedure. Do you think that is possible? Do you have any suggestion on how to do it?

Thank you again, Sebastian

— You are receiving this because you commented. Reply to this email directly, view it on GitHub

<https://github.com/harrelfe/rms/issues/26#issuecomment-229942984 , or mute the thread

<

https://github.com/notifications/unsubscribe/ABGO2nth7gJtvl6G8kgAp6jbNtUtfl-gks5qRRJXgaJpZM4JCW81

.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub

https://github.com/harrelfe/rms/issues/26#issuecomment-229952667, or mute the thread <

https://github.com/notifications/unsubscribe/AFJu9vkIeunp1aqUBPEex_EBTYULC5FKks5qRRvWgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-229992802, or mute the thread

<

https://github.com/notifications/unsubscribe/ABGO2kGBdYviFLDokCNARXSbdp1LHhHDks5qRULVgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230009008, or mute the thread <

https://github.com/notifications/unsubscribe/AFJu9tLYLeuyG2JRx3JY9H1fH1bPd2woks5qRVPLgaJpZM4JCW81

.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230010835, or mute the thread

https://github.com/notifications/unsubscribe/ABGO2jh56KD4OIMkwvn5bytAGd1RP01kks5qRVXCgaJpZM4JCW81 .

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/harrelfe/rms/issues/26#issuecomment-230038111, or mute the thread https://github.com/notifications/unsubscribe/AFJu9mv0M6vehGC0EE8qsPEOoLLCR9Rnks5qRXWAgaJpZM4JCW81.