[R-sig-ME] Predicting values from an MCMCglmm model

[Ashby, Matthew]

Dear all,

I am fitting an MCMCglmm with Poisson link to investigate relationships between schools and violence in the surrounding neighbourhood. The dependent variable is counts of assaults in neighbourhoods, with repeated measures across four different time periods. The model uses data from nine cities. The model is fit by:

aslt_h1 <- MCMCglmm(
  fixed = count ~ count_lag_sc,
  random = ~ us(1 + pop_total_sc + index_disadvantage + index_mobility + index_ethnic + perc_teen_sc + private_elementary_school + public_elementary_school + private_secondary_school +  public_secondary_school):city_period,
  data = data_h1,
  family = "poisson",
  prior = list(R = list(V = diag(1), nu = 0.002), G = list(G1 = list(V = diag(10), n = 10))),
  thin = 20, burnin = 5000, nitt = 50000, pr = TRUE
)

I would like to generate predicted values for the mean neighbourhood at each time period in each city. I've constructed a data frame of new data, but I get the error message

Error in dimnames(x) <- dn :
  length of 'dimnames' [2] not equal to array extent

By trial and error I think the problem is caused by my newdata data frame not being in the format expected by predict.MCMCglmm(), but having read several tutorials online and section 2.5 of the Course Notes, I can't see where I'm going wrong. I'm creating the newdata DF by:

newdata <- expand.grid(
  count = 0,
  # since all the control variables were centered and scaled, they can be set to zero here to represent block groups with the mean values
  count_lag_sc = 0,
  pop_total_sc = 0,
  perc_teen_sc = 0,
  index_disadvantage = 0,
  index_mobility = 0,
  index_ethnic = 0,
  # since we know there is only a meaningful effect for public high/middle schools, we can set all the other school indicators to FALSE
  private_elementary_school = FALSE,
  public_elementary_school = FALSE,
  private_secondary_school = FALSE,
  public_secondary_school = c(FALSE, TRUE),
  city_period = c(
    "Chicago - non-school day non-school time",
    "Chicago - non-school day school time",
    "Chicago - school day non-school time",
    "Chicago - school day school time",
    "Detroit - non-school day non-school time",
    [… levels for other cities snipped for space …]
  )
)

And then trying to predict values by calling:

predictions <- predict(aslt_h1, newdata = newdata, marginal = NULL, type = "response", posterior = "mean", verbose = TRUE)

Does the format of the newdata DF look acceptable? Is the problem likely to be there or in the arguments I've set in the predict() call?

Thank you in advance for any tips or suggestions.

Yours,


Matt Ashby
Nottingham Trent University


DISCLAIMER: This email is intended solely for the addressee. It may contain private and confidential information. If you are not the intended addressee, please take no action based on it nor show a copy to anyone. In this case, please reply to this email to highlight the error. Opinions and information in this email that do not relate to the official business of Nottingham Trent University shall be understood as neither given nor endorsed by the University. Nottingham Trent University has taken steps to ensure that this email and any attachments are virus-free, but we do advise that the recipient should check that the email and its attachments are actually virus free. This is in keeping with good computing practice.