Chapter 11 - God Spiked the Integers

This chapter described some of the most common generalized linear models, those used to model counts. It is important to never convert counts to proportions before analysis, because doing so destroys information about sample size. A fundamental difficulty with these models is that parameters are on a different scale, typically log-odds (for binomial) or log-rate (for Poisson), than the outcome variable they describe. Therefore computing implied predictions is even more important than before.

Place each answer inside the code chunk (grey box). The code chunks should contain a text response or a code that completes/answers the question or activity requested. Make sure to include plots if the question requests them. Problems are labeled Easy (E), Medium (M), and Hard(H).

Finally, upon completion, name your final output .html file as: YourName_ANLY505-Year-Semester.html and publish the assignment to your R Pubs account and submit the link to Canvas. Each question is worth 5 points.

Questions

11E1. If an event has probability 0.35, what are the log-odds of this event?

options(digits=5)
P <- 0.35
log_odds <- log(P/(1-P))

11E2. If an event has log-odds 3.2, what is the probability of this event?

# log_odds = log(P/(1-P))
# exp(log_odds) = P/(1-P)
# P = exp(log_odds)/(1+exp(log_odds))
P <- exp(3.2)/(1+exp(3.2))
P
## [1] 0.96083

11E3. Suppose that a coefficient in a logistic regression has value 1.7. What does this imply about the proportional change in odds of the outcome?

coefficient <- exp(1.7)
coefficient
## [1] 5.4739
# This implies that we are 5.5 times more likely to reach the given outcome.

11E4. Why do Poisson regressions sometimes require the use of an offset? Provide an example.

# The offset, also known as exposure, is used as a correction in Poisson regressions to provide a uniform setting for events happening is different scales of time/distance. An example could be predicting the number of car crashes in a month vs in a year.

11M1. As explained in the chapter, binomial data can be organized in aggregated and disaggregated forms, without any impact on inference. But the likelihood of the data does change when the data are converted between the two formats. Can you explain why?

# The reason behind the difference in the likelihood between the two formats lies in the presence of an extra log odd factor in the aggregated form. As the conversion happens to the diaggregated form, the multiplier c(n,m) turns to a constant, thus changing the likelihood outcome.

11M2. If a coefficient in a Poisson regression has value 1.7, what does this imply about the change in the outcome?

coeff <- 1.7
exp(coeff)
## [1] 5.4739
# Every time the predictor variable changes by 1 unit, the lambda parameter of the Poisson regression will increase 5.47 times.

11M3. Explain why the logit link is appropriate for a binomial generalized linear model.

# Given the fact that the outcome generated by a binominal generalized model is consistently binary, the logit link should be appropriate for it since it maps the probability of a given variable as a parameter with a value limited between 0 and 1.

11M4. Explain why the log link is appropriate for a Poisson generalized linear model.

# As motioned in the answer above, the log link's predictor value is positive, between 0 and 1, which satisfies the requirements of a Poisson generalized linear model that produces a positive outcome at all times.

11M5. What would it imply to use a logit link for the mean of a Poisson generalized linear model? Can you think of a real research problem for which this would make sense?

# It implies that the outcome likelihood will always be positive. I unfortunately wasn't able to come up with a research problem to illustrate.

11M6. State the constraints for which the binomial and Poisson distributions have maximum entropy. Are the constraints different at all for binomial and Poisson? Why or why not?

# The constraints for which the binomial and Poisson distributions have maximum entropy are, respectively: 
# - Events are discrete.
# - Expected value is constant.
# Given the fact that the distributions are different, the constraints are different as well. for the Poisson distribution, the variance is constant and must be equal to expected values.

11M7. Use quap to construct a quadratic approximate posterior distribution for the chimpanzee model that includes a unique intercept for each actor, m11.4 (page 330). Plot and compare the quadratic approximation to the posterior distribution produced instead from MCMC. Can you explain both the differences and the similarities between the approximate and the MCMC distributions? Relax the prior on the actor intercepts to Normal(0,10). Re-estimate the posterior using both ulam and quap. Plot and compare the posterior distributions. Do the differences increase or decrease? Why?

data("chimpanzees")
data <- chimpanzees
data$recipient <- NULL

m <- map(alist(
  pulled_left ~ dbinom( 1 , p ) ,
  logit(p) <- a[actor] + (bp + bpC*condition)*prosoc_left ,
  a[actor] ~ dnorm(0,10),
  bp ~ dnorm(0,10),
  bpC ~ dnorm(0,10)
) ,
data=data)
pairs(m)

# We can observe that the posterior standard deviation is almost similar to the posterior mean. 
# We can also make the remark that the MCMC model's posterior standard deviation is a little higher.

11M8. Revisit the data(Kline) islands example. This time drop Hawaii from the sample and refit the models. What changes do you observe?

data(Kline)
data2 <- Kline
data2$P <- scale( log(data2$population) )
data2$contact_id <- ifelse( data2$contact=="high" , 2 , 1 )
data2
##       culture population contact total_tools mean_TU          P contact_id
## 1    Malekula       1100     low          13     3.2 -1.2914733          1
## 2     Tikopia       1500     low          22     4.7 -1.0885508          1
## 3  Santa Cruz       3600     low          24     4.0 -0.5157649          1
## 4         Yap       4791    high          43     5.0 -0.3287734          2
## 5    Lau Fiji       7400    high          33     5.0 -0.0443390          2
## 6   Trobriand       8000    high          19     4.0  0.0066683          2
## 7       Chuuk       9200    high          40     3.8  0.0981092          2
## 8       Manus      13000     low          28     6.6  0.3243176          1
## 9       Tonga      17500    high          55     5.4  0.5187979          2
## 10     Hawaii     275000     low          71     6.6  2.3210083          1

11H1. Use WAIC or PSIS to compare the chimpanzee model that includes a unique intercept for each actor, m11.4 (page 330), to the simpler models fit in the same section. Interpret the results.

data3 <- chimpanzees

m11.1 <- map(
  alist(
    pulled_left ~ dbinom(1, p),
    logit(p) <- a ,
    a ~ dnorm(0,10)
  ),
  data=data3 )


m11.2.2 <- map(
  alist(
    pulled_left ~ dbinom(1, p) ,
    logit(p) <- a + bp*prosoc_left ,
    a ~ dnorm(0,10) ,
    bp ~ dnorm(0,10)
  ),
  data=data3 )

m11.3 <- map(
  alist(
    pulled_left ~ dbinom(1, p) ,
    logit(p) <- a + (bp + bpC*condition)*prosoc_left ,
    a ~ dnorm(0,10) ,
    bp ~ dnorm(0,10) ,
    bpC ~ dnorm(0,10)
  ), data=data3 )

m11.4 <- map(
  alist(
    pulled_left ~ dbinom(1, p),
    logit(p) <- a[actor] + (bp + bpC*condition)*prosoc_left,
    a[actor] ~ dnorm(0, 10),
    bp ~ dnorm(0, 10),
    bpC ~ dnorm(0, 10)
  ),
  data = data3)

compare(m11.1,m11.2.2,m11.3,m11.4)
##           WAIC      SE  dWAIC    dSE   pWAIC     weight
## m11.4   542.50 19.0047   0.00     NA 12.1343 1.0000e+00
## m11.2.2 680.42  9.3514 137.92 18.349  1.9626 1.1250e-30
## m11.3   682.32  9.3916 139.81 18.277  2.9880 4.3620e-31
## m11.1   687.97  7.0996 145.47 19.150  1.0161 2.5793e-32
# m11.4 appears to be the best model out of all the above.