## Model: Linear regression
## Fit: lm
## Data: 1920s cars data
data("cars", package = "datasets")
## Stopping distance (ft) explained by speed (mph)
reg <- lm(dist ~ speed, data = cars)
## Extract fitted normal distributions (in-sample, with constant variance)
pd <- prodist(reg)
head(pd)
## Extract log-likelihood from model object
logLik(reg)
## Replicate log-likelihood via distributions object
sum(log_pdf(pd, cars$dist))
log_likelihood(pd, cars$dist)
## Compute corresponding medians and 90% interval
qd <- quantile(pd, c(0.05, 0.5, 0.95))
head(qd)
## Visualize observations with predicted quantiles
plot(dist ~ speed, data = cars)
matplot(cars$speed, qd, add = TRUE, type = "l", col = 2, lty = 1)
## Sigma estimated by maximum-likelihood estimate (default, used in logLik)
## vs. least-squares estimate (used in summary)
nd <- data.frame(speed = 50)
prodist(reg, newdata = nd, sigma = "ML")
prodist(reg, newdata = nd, sigma = "OLS")
summary(reg)$sigma
## Model: Poisson generalized linear model
## Fit: glm
## Data: FIFA 2018 World Cup data
data("FIFA2018", package = "distributions3")
## Number of goals per team explained by ability differences
poisreg <- glm(goals ~ difference, data = FIFA2018, family = poisson)
summary(poisreg)
## Interpretation: When the ratio of abilities increases by 1 percent,
## the expected number of goals increases by around 0.4 percent
## Predict fitted Poisson distributions for teams with equal ability (out-of-sample)
nd <- data.frame(difference = 0)
prodist(poisreg, newdata = nd)
## Extract fitted Poisson distributions (in-sample)
pd <- prodist(poisreg)
head(pd)
## Extract log-likelihood from model object
logLik(poisreg)
## Replicate log-likelihood via distributions object
sum(log_pdf(pd, FIFA2018$goals))
log_likelihood(pd, FIFA2018$goals)
## Model: Autoregressive integrated moving average model
## Fit: arima
## Data: Quarterly approval ratings of U.S. presidents (1945-1974)
data("presidents", package = "datasets")
## ARMA(2,1) model
arma21 <- arima(presidents, order = c(2, 0, 1))
## Extract predicted normal distributions for next two years
p <- prodist(arma21, n.ahead = 8)
p
## Compute median (= mean) forecast along with 80% and 95% interval
quantile(p, c(0.5, 0.1, 0.9, 0.025, 0.975))
Run the code above in your browser using DataLab