I've developed a glmer model to predict the presence/absence of a certain behavior based on multiple environmental predictors (temperature, windspeed, precipitation). 14 animals were observed continuously throughout the daylight hours over multiple days, and data were summarized over two-hour periods. The response is the proportion of the two-hour period that the individual spent performing the behavior, and the predictors are the average temperature, windspeed, and precipitation over that two-hour period. Because the response variable is a proportion, I've included a weighting factor in the model. I'm most interested in the relationship between temperature and the response behavior. The model is specified as:

Model <- glmer(Response~Temp+Precip+Windspeed+(Temp|Individual),
   data=data, family=binomial, weights=weightfactor)

I'd like to calculate the temperature at which, with other predictors (precipitation and wind speed) held at their mean values, the probability of the response = 0.5, and I'd also like to calculate 95% confidence intervals for that estimate. I found these threads, which describe how to adapt the dose.p function to work with glmer models. I've used those as a basis for my code:

  dose.p.glmm <-  function(Model, cf = 1:2, p = 0.5) {
  f <- family(Model)
  eta <- f$linkfun(p)
  b <- fixef(Model)[cf]
  x.p <- (eta - b[1L])/b[2L]
  names(x.p) <- paste("p = ", format(p), ":", sep = "")
  pd <- -cbind(1, x.p)/b[2L]
  SE <- sqrt(((pd %*% vcov(Model)[cf, cf]) * pd) %*% c(1, 1))
  res <- structure(x.p, SE = matrix(SE), p = p)
  class(res) <- "glm.dose"
  res}
  dose.p.glmm(Model, cf=1:2, p=0.5)

However, this code was written for glmer models with a single predictor variable. Is there a way to adapt my code to specify that in calculating this, I want to hold precip and windspeed at their mean observed values and then find the value of temp at which the probability of the response=0.5? And if not, is there another way to calculate that value and its 95% confidence intervals?