-2

I have the following model : 

model_ <- glm( response ~ var_1 + var_2, family = "binomial" )

which gives me the following results :

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) -3.07418    0.05484 -56.053   <2e-16 ***
var_1        0.19238    0.13547   1.420    0.156    
var_2        2.07090    0.23579   8.783   <2e-16 ***
---

The only way I figured out how to plot with the base graphic functions the two probability curves was by running two different models :

model_1 <- glm( response ~ var_1, family = "binomial" )
model_2 <- glm( response ~ var_2, family = "binomial" )

and create two predicted curves :

curve( predict( mod_1, data.frame( var_1=x ), type="resp" ))
curve( predict( mod_2, data.frame( var_2=x ), type="resp" ))

Which once plotted look like the following graphs :

plots_logreg

But of course the coefficients are then not the same than with my initial model_. 

          coef(model_)    coef(model_1)    coef(model_2)
var_1         0.192            0.165           ---
var_2         2.071***          ---           2.065***

I would like to plot the two probability curves corresponding to each variable from the multivariate model model_and not from the two submodels model_1and model_2.

  • Is there a more or less simple way to show those two curves directly from model_ without having to run other submodels which give slightly different coefficients ?
  • How would it work with 4 probability curves corresponding to the variables var_1, var_2, var_3and var_4 with a formula such as response ~ var_1 + var_2 + var_3 + var_4 ?

There are some examples on the net, but either using automatic functions from ggplot (which i don't want to use) or running univariate models.

Adrx
  • 31
  • 1
  • 9
  • 2
    can you `dput` your sample data? – Sandipan Dey Mar 13 '17 at 15:55
  • My data are 5 vectors of 8800 observations each. I can't make them to ASCII here. – Adrx Mar 13 '17 at 16:00
  • 1
    When asking for help, you should include a [reproducible example](http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) with sample input and code we can run. This will make it easier to help you. Also be clear on exactly what you want this plot to look like. Plots only have one x and y axis at a time. If you are looking for visualizations recommendations rather than asking a specific programming question, you might want to ask at [stats.se] instead. – MrFlick Mar 13 '17 at 16:04
  • There are packages that can help you, such as the `effects` package, which can easily give you `y ~ x` predictions on the mean of covariates. – Axeman Mar 13 '17 at 16:11
  • The plot i have uploaded runs with two submodels (`model_1`and `model_2`), one for each variable `var_1`and `var_2`. They give different coefficients as with only one multivariate model (`model_`). I would like to correct these curves plot them with the coefficients given by the multivariate model, and not by the two submodels. – Adrx Mar 13 '17 at 16:16

1 Answers1

0

I came out with this last night : 

# intercept & beta coefficients  
intercept <- as.numeric( coef( model )[1] )
beta_var_1 <- as.numeric( coef( model )[2] )
beta_var_2 <- as.numeric( coef( model )[3] )
...

# range of variables
range <- seq( 0, 10, .01 )

# log odds
logodds_var_1 <- intercept + beta_var_1 * range
logodds_var_2 <- intercept + beta_var_2 * range
...

# probabilities 
prob_var_1 <- exp(logodds_var_1) / (1+exp(logodds_var_1))
prob_var_2 <- exp(logodds_var_2) / (1+exp(logodds_var_2))
...

# curves
plot( x=range, y=prob_var_1, type="l" )
plot( x=range, y=prob_var_2, type="l" )
...

This produces two (or more) curves which look more coherent to what i am looking for. I hope i am not misunderstanding the principles of such a logit model.

Adrx
  • 31
  • 1
  • 9