Skip to content

Instantly share code, notes, and snippets.

@ito4303
Created December 22, 2022 06:50
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save ito4303/4509fafa5390d167c050df5733c9f612 to your computer and use it in GitHub Desktop.
Save ito4303/4509fafa5390d167c050df5733c9f612 to your computer and use it in GitHub Desktop.
Likelihood ratio test in R (using base functions or lrtest function)
---
title: "R Notebook"
output: html_notebook
---
## An example of likelihood ratio tests in R
```{r setup}
library(lmtest)
library(palmerpenguins)
library(ggplot2)
```
## View data
Adelie penguin data was used.
Omitted the missing values assuming MCAR.
```{r}
adelie <- penguins[penguins$species == "Adelie" &
!is.na(penguins$bill_depth_mm) &
!is.na(penguins$bill_length_mm), ]
ggplot(adelie) +
geom_point(aes(x = bill_depth_mm, y = bill_length_mm, colour = island))
ggsave("plot01.png", width = 640, height = 480, units = "px", dpi = 144)
```
### Models
Model1: bill_length_mm ~ bill_depth_mm
Model2: bill_length_mm ~ bill_depth_mm + island
```{r}
model1 <- lm(bill_length_mm ~ bill_depth_mm, data = adelie)
model2 <- lm(bill_length_mm ~ bill_depth_mm + island, data = adelie)
```
### Likelihood ratio test
using the base functions
```{r}
(logLik1 <- logLik(model1))
(logLik2 <- logLik(model2))
df1 <- attr(logLik1, "df")
df2 <- attr(logLik2, "df")
d2logLik <- 2 * (logLik2 - logLik1)
attributes(d2logLik) <- NULL
d2logLik
1 - pchisq(d2logLik, df = df2 - df1)
```
using the lrtest function
```{r}
lrtest(model1, model2)
```
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment