Search code examples
rglmbinomial-coefficients

Plotting binomial glm with interactions in numeric variables


I want to know if is possible to plotting binomial glm with interactions in numeric variables. In my case:

##Data set artificial
set.seed(20)
d <- data.frame(
    mating=sample(0:1, 200, replace=T),
    behv = scale(rpois(200,10)),
    condition = scale(rnorm(200,5))
) 

#Binomial GLM ajusted
model<-glm(mating ~ behv + condition, data=d, family=binomial)
summary(model)

In a situation where behv and condition are significant in the model

#Plotting first for behv
x<-d$behv ###Take behv values
x2<-rep(mean(d$condition),length(d_p[,1])) ##Fixed mean condition

# Points
plot(d$mating~d$behv)

#Curve
curve(exp(model$coefficients[1]+model$coefficients[2]*x+model$coefficients[3]*x2)
/(1+exp(model$coefficients[1]+model$coefficients[2]*x+model$coefficients[3]*x2)))

But doesn't work!! There is another correct approach?

Thanks


Solution

  • It seems like your desired output is a plot of the conditional means (or best-fit line). You can do this by computing predicted values with the predict function.

    I'm going to change your example a bit, to get a nicer looking result.

    d$mating <- ifelse(d$behv > 0, rbinom(200, 1, .8), rbinom(200, 1, .2))
    model <- glm(mating ~ behv + condition, data = d, family = binomial)
    summary(model)
    

    Now, we make a newdata dataframe with your desired values:

    newdata <- d
    newdata$condition <- mean(newdata$condition)
    newdata$yhat <- predict(model, newdata, type = "response")
    

    Finally, we sort newdata by the x-axis variable (if not, we'll get lines that zig-zag all over the plot), and then plot:

    newdata <- newdata[order(newdata$behv), ]
    plot(newdata$mating ~ newdata$behv)
    lines(x = newdata$behv, y = newdata$yhat)
    

    Output:

    enter image description here