Interferensi dalam Regresi Linier

1 Interferensi dalam Regresi Linier

Setelah membaca bab ini, anda akan dapat:

• Memahami distribusi estimasi regresi.
• Membuat interval untuk parameter regresi, responden rata-rata, dan prediksi.
• Menguji taraf signifikansi suatu regresi.

Dalam bab sebelumnya, anda sudah menegaskan model dari regresi linear sederhana, \[Y_i = \beta_0 + \beta_1 x_i + \epsilon_i\] dimana \(\epsilon_i \sim N(0, \sigma^2)\). Lalu gunakan observasi \((x_i,y_i)\), dengan \(i=1,2,...,n\) untuk mendapatkan nilai \(\beta_0\) dan \(\beta_1\) yang diminimalkan \[f(\beta_0, \beta_1) = \sum_{i = 1}^{n}(y_i - (\beta_0 + \beta_1 x_i))^2.\] Nilainya dapat juga disebut \(\hat{\beta_0}\) dan \(\hat{\beta_1}\), dapat dibentuk juga menjadi \[\begin{aligned} \hat{\beta}_1 &= \frac{S_{xy}}{S_{xx}} = \frac{\sum_{i = 1}^{n}(x_i - \bar{x})(y_i - \bar{y})}{\sum_{i = 1}^{n}(x_i - \bar{x})^2}\\ \hat{\beta}_0 &= \bar{y} - \hat{\beta}_1 \bar{x}. \end{aligned}\] \(\sigma^2\) juga di estimasi menggunakan \(s_e^2\). Dalam kata lain, \(s_e\) merupakan estimasi dari \(\sigma\), dimana \[s_e = \text{RSE} = \sqrt{\frac{1}{n - 2}\sum_{i = 1}^n e_i^2}\] dapat disebut RSE (Residual Standard Error).

Ketika diaplikasikan dalam data cars, berikut hasil yang diperoleh:

stop_dist_model = lm(dist ~ speed, data = cars)
summary(stop_dist_model)

## 
## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -29.069  -9.525  -2.272   9.215  43.201 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
## speed         3.9324     0.4155   9.464 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.38 on 48 degrees of freedom
## Multiple R-squared:  0.6511, Adjusted R-squared:  0.6438 
## F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

Bagian sebelumnya, kita hanya membahas Estimasi, RSE, dan R-kuadrat berganda. Kali ini akan membahas mengenai koefisien dan F-statistik.

plot(dist ~ speed, data = cars,
     xlab = "Speed (in Miles Per Hour)",
     ylab = "Stopping Distance (in Feet)",
     main = "Stopping Distance vs Speed",
     pch  = 20,
     cex  = 2,
     col  = "grey")
abline(stop_dist_model, lwd = 5, col = "darkorange")

Sebelum mulai, berikut catatan untuk persamaan \(S_{xy}\) yang tidak diberikan pada bagian sebelumnya, \[S_{xy}= \sum_{i = 1}^{n}(x_i - \bar{x})(y_i - \bar{y}) = \sum_{i = 1}^{n}(x_i - \bar{x}) y_i.\] Persamaan ini mungkin tidak terduga. (Coba untuk buktikan.) Namun, akan berguna untuk mengilustrasikan konsep pada bagian ini.

Perlu diingat bahwa \(\hat{\beta_1}\) adalah sampel statistik ketika dikalkulasikan dengan data yang diperlihatkan diatas seperti \(\hat{\beta_0}\)

Namun, dalam bagian ini wajar untuk menggunakan \(\hat{\beta_1}\) dan \(\hat{\beta_0}\) sebagai variabel acak, yang berarti nilai \(Y_i\) belum diperoleh. Dalam kasus ini, notasi yang akan digunakan sedikit berbeda, mensubstitusi \(Y_i\) dan \(y_i\). \[\begin{aligned} \hat{\beta}_1 &= \frac{\sum_{i = 1}^{n}(x_i - \bar{x}) Y_i}{\sum_{i = 1}^{n}(x_i - \bar{x})^2} \\ \hat{\beta}_0 &= \bar{Y} - \hat{\beta}_1 \bar{x} \end{aligned}\] Pada bagian sebelumnya, dinyatakan bahwa parameter \(\beta_0\) dan \(\beta_1\) untuk model yang tidak diketahui dapat diperoleh dengan meminimalkan error. Sekarang membahas teorema Gauss-Markov yang melanjutkan permasalahan tersebut lebih dalam lagi, memperlihatkan beberapa estimasi yang berupa estimasi “terbaik” dalam sudut pandang tertentu.

1.1 Teorema Gauss-Markov

Teorema Gauss-Markov mengatakan bahwa saat melakukan estimasi parameter untuk regresi linear sederhana \(\beta_0\),\(\beta_1\), \(\beta_0\) dan \(\beta_1\) yang diturunkan berupa best linear unbiased estimates (BLUE). (Kondisi sebenarnya untuk teorema Gauss-Markov lebih fleksibel dibandingan dengan model SLR)

Sekarang membahas linier, unbiased, dan terbaik yang bersangkutan dengan estimasi ini.

1.1.1 Linier

Mengingat kembali aturan SLR yaitu nilai \(x_i\) dikatakan tetap dengan nilai yang diketahui. Maka estimasi linier dapat dituliskan sebagai kombinasi linier dari \(Y_i\). Dalam kasus dari \(\hat{\beta_1}\) dapat dilihat \[\hat{\beta}_1 = \frac{\sum_{i = 1}^{n}(x_i - \bar{x}) Y_i}{\sum_{i = 1}^{n}(x_i - \bar{x})^2} = \sum_{i = 1}^n k_i Y_i = k_1 Y_1 + k_2 Y_2 + \cdots k_n Y_n\] dimana \(k_i = \displaystyle\frac{(x_i - \bar{x})}{\sum_{i = 1}^{n}(x_i - \bar{x})^2}\)

Dalam bentuk yang mirip, dapat dilihat bahwa \(\hat{\beta_0}\) dapat dituliskan sebagai kombinasi linier dari \(Y_i\). Maka \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) adalah estimator linier.

1.1.2 Unbiased

Sekarang telah diketahui bahwa estimasinya linier, seberapa baguskah estimasinya? Salah satu perhitungan “goodness” suatu estimate adalah bias. Secara spesifik, lebih baik jika estimasi tersebut tidak bias, yang artinya nilai ekspektasinya adalah parameter yang diestimasikan.

Dalam kasus estimasi regresi, diperoleh: \[\begin{aligned} \text{E}[\hat{\beta}_0] &= \beta_0 \\ \text{E}[\hat{\beta}_1] &= \beta_1. \end{aligned}\]

Dengan ini dapat dikatakan bahwa saat kondisi SLR terpenuhi, maka estimasinya rata-rata benar. Namun, seperti pada bagian sebelumnya saat mensimulasikan model SLR itu tidak menunjukkan jika estimasi individunya benar. Jika proses tersebut diulang sampai dengan tak terhingga pengulangan, maka estimasinya rata-rata benar.

1.1.3 Terbaik

Sekarang, jika estimasi linier dan unbiased dibatasi, bagaimana cara menjelaskan estimasi terbaik? Dengan variansi minimum.

Pertama, perlu dicatat bahwa sangat mudah untuk membuat estimasi untuk \(\beta_1\) yang memiliki variansi sangat rendah tetapi bukan unbiased. Contohnya: \[\hat{\theta}_{BAD} = 5.\] Karena \(\hat{\theta}_{BAD}\) adalah nilai konstan \[\text{Var}[\hat{\theta}_{BAD}] = 0.\] Namun sejak, \[\text{E}[\hat{\theta}_{BAD}] = 5\] \(\hat{\theta}_{BAD}\) dapat dikatakan estimator bias kecuali \(\beta_1=5\), yang berarti tidak diketahui sebelumnya. Untuk hal ini, estimasinya buruk (kecuali \(\beta_1=5\)) walaupun variansi yang terkecil. Ini adalah alasan mengapa estimasi unbiased tersebut dibatasi. Apa bagusnya suatu estimasi jika mengestimasi jumlah yang tidak tepat?

Maka, pertanyaannya adalah, apa variansi dari \(\hat{\beta_0}\) dan \(\hat{\beta_1}\)?

\[\begin{aligned} \text{Var}[\hat{\beta}_0] &= \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}\right) \\ \text{Var}[\hat{\beta}_1] &= \frac{\sigma^2}{S_{xx}}. \end{aligned}\]

Hal ini mengukur variabilitas dari suatu estimasi dikarenakan kemungkinan yang acak saat sampling. Apakah ini yang “terbaik”? Apakah variansinya sudah yang paling kecil yang dapat diperoleh? Anda hanya perlu mempercayai apa yang kami katakan karena menunjukkan bahwa hal itu benar berada di luar cakupan kursus ini.

1.2 Distribusi sampling

Sekarang estimasi \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) sebagai variabel acak, maka dapat didiskusikan distribusi samplingnya yang berarti distribusi saat statistik dianggap variabel acak.

Karena \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) berupa kombinasi linier dari \(Y_i\) dan tiap \(Y_i\) distribusi normal, maka \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) juga mengikuti distribusi normal.

Lalu dengan menggabungkan semua menjadi satu menghasilkan distribusi \(\hat{\beta_0}\) dan \(\hat{\beta_1}\).

Untuk \(\hat{\beta_1}\), \[\hat{\beta}_1 = \frac{S_{xy}}{S_{xx}} = \frac{\sum_{i = 1}^{n}(x_i - \bar{x}) Y_i}{\sum_{i = 1}^{n}(x_i - \bar{x})^2} \sim N\left( \beta_1, \ \frac{\sigma^2}{\sum_{i = 1}^{n}(x_i - \bar{x})^2} \right).\]

Jika disederhanakan, \[\hat{\beta}_1 \sim N\left( \beta_1, \frac{\sigma^2}{S_{xx}} \right).\] dan untuk \(\hat{\beta_0}\), \[\hat{\beta}_0 = \bar{Y} - \hat{\beta}_1 \bar{x} \sim N\left( \beta_0, \ \frac{\sigma^2 \sum_{i = 1}^{n}x_i^2}{n \sum_{i = 1}^{n}(x_i - \bar{x})^2} \right).\]

Jika disederhanakan, \[\hat{\beta}_0 \sim N\left( \beta_0, \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}\right) \right)\]

Saat ini belum ada bukti untuk hasilnya. Hasil tersebut akan dicocokan dengan menggunakan simulasi dan bukan melakukan turunan dari distribusi sampling.

Catatan untuk pembaca: Turunan dan bukti dapat ditambahkan dalam lampiran suatu waktu nanti. Anda juga bisa menemukan hasil yang mendekati pada buku regresi standar linier. Di UIUC, hasil tersebut kemungkinan dipresentasikan pada STAT 424 dan STAT 425. Namun, karena anda tidak akan mempertanyakan penurunan, dalam kursus ini tidak akan dibahas.

1.2.1 Simulasi Distribusi Sampling

Untuk membuktikan hasil diatas, simulasikan sample ukuran \(n=100\) dari model \[Y_i = \beta_0 + \beta_1 x_i + \epsilon_i\] dimana \(\epsilon_i \sim N(0, \sigma^2).\) Dalam kasus ini parameternya:

• \(\beta_0=3\)
• \(\beta_1=6\)
• \(\sigma^2=4\)

Dari pernyataan diatas, menghasilkan \[\hat{\beta}_1 \sim N\left( \beta_1, \frac{\sigma^2}{S_{xx}} \right)\] dan \[\hat{\beta}_0 \sim N\left( \beta_0, \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}\right) \right).\] Pertama tentukan dahulu nilai \(x\) untuk simulasi, karena \(x\) dalam SLR juga dianggap kuantitas. Penentuan \(x\) tidak diatur. Disini juga mempersiapkan seed untuk pengacakan, dan kalkulasi \(S_{xx}\) yang diperlukan nantinya.

set.seed(42)
sample_size = 100 # this is n
x = seq(-1, 1, length = sample_size)
Sxx = sum((x - mean(x)) ^ 2)

Memperbaiki nilai parameter

beta_0 = 3
beta_1 = 6
sigma  = 2

Dengan informasi diatas, distribusi samplingnya:

(var_beta_1_hat = sigma ^ 2 / Sxx)

## [1] 0.1176238

(var_beta_0_hat = sigma ^ 2 * (1 / sample_size + mean(x) ^ 2 / Sxx))

## [1] 0.04

\[\hat{\beta}_1 \sim N( 6, 0.1176238)\] dan \[\hat{\beta}_0 \sim N( 3, 0.04).\] Hasilnya, \[\begin{aligned} \text{E}[\hat{\beta}_1] &= 6 \\ \text{Var}[\hat{\beta}_1] &= 0.1176238 \end{aligned}\] dan \[\begin{aligned} \text{E}[\hat{\beta}_0] &= 3 \\ \text{Var}[\hat{\beta}_0] &= 0.04. \end{aligned}\]

Sekarang simulasikan data dari model ini sebanyak 10.000 kali. Catatan: Hal ini bukan cara R untuk melakukan simulasi. Simulasi dengan cara ini dilakukan untuk menambah kejelasan. Contoh: Fungsi sim_slr() dapat digunakan pada bagian sebelumnya. Variabel tersebut disimpan pada global environment daripada membuat data frame untuk tiap simulasi.

num_samples = 10000
beta_0_hats = rep(0, num_samples)
beta_1_hats = rep(0, num_samples)

for (i in 1:num_samples) {
  eps = rnorm(sample_size, mean = 0, sd = sigma)
  y   = beta_0 + beta_1 * x + eps
  
  sim_model = lm(y ~ x)
  
  beta_0_hats[i] = coef(sim_model)[1]
  beta_1_hats[i] = coef(sim_model)[2]
}

Tiap simulasi data dilakukan menghasilkan koefisien estimasi. Variabel beta_0_hats dan beta_1_hats sekarang menyimpan 10.000 nilai simulasi dari \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) secara berurutan.

Pertama membuktikan distribusi dari \(\hat{\beta_1}\).

mean(beta_1_hats) # empirical mean

## [1] 6.001998

beta_1            # true mean

## [1] 6

var(beta_1_hats)  # empirical variance

## [1] 0.11899

var_beta_1_hat    # true variance

## [1] 0.1176238

Nilai empiris rata-rata dan variansi sangat mirip. Sudah dibuktikan distribusinya itu normal. beta_1_hats divisualisasikan dalam bentuk histogram, dan menambahkan kurva distribusi \(\hat{\beta_1}\). Untuk menyesuaikan skala dengan kurva menggunakan prob = TRUE.

# note need to use prob = TRUE
hist(beta_1_hats, prob = TRUE, breaks = 20, 
     xlab = expression(hat(beta)[1]), main = "", border = "dodgerblue")
curve(dnorm(x, mean = beta_1, sd = sqrt(var_beta_1_hat)), 
      col = "darkorange", add = TRUE, lwd = 3)

Lalu ulang proses \(\hat{\beta_0}\).

mean(beta_0_hats) # empirical mean

## [1] 3.001147

beta_0            # true mean

## [1] 3

var(beta_0_hats)  # empirical variance

## [1] 0.04017924

var_beta_0_hat    # true variance

## [1] 0.04

hist(beta_0_hats, prob = TRUE, breaks = 25, 
     xlab = expression(hat(beta)[0]), main = "", border = "dodgerblue")
curve(dnorm(x, mean = beta_0, sd = sqrt(var_beta_0_hat)),
      col = "darkorange", add = TRUE, lwd = 3)

Pembelajaran simulasi ini hanya mensimulasikan jumlah sampel yang terbatas. Untuk membuktikan hasil distribusi, diperlukan sampel yang tidak terhingga. Namun dalam plot berikut memperlihatkannya dengan cukup jelas jika simulasi itu berlanjut.

par(mar = c(5, 5, 1, 1)) # adjusted plot margins, otherwise the "hat" does not display
plot(cumsum(beta_1_hats) / (1:length(beta_1_hats)), type = "l", ylim = c(5.95, 6.05),
     xlab = "Number of Simulations",
     ylab = expression("Empirical Mean of " ~ hat(beta)[1]),
     col  = "dodgerblue")
abline(h = 6, col = "darkorange", lwd = 2)

par(mar = c(5, 5, 1, 1)) # adjusted plot margins, otherwise the "hat" does not display
plot(cumsum(beta_0_hats) / (1:length(beta_0_hats)), type = "l", ylim = c(2.95, 3.05),
     xlab = "Number of Simulations",
     ylab = expression("Empirical Mean of " ~ hat(beta)[0]),
     col  = "dodgerblue")
abline(h = 3, col = "darkorange", lwd = 2)

1.3 Kesalahan Baku (Standard Error)

Sekarang ada 2 hasil distribusi yang dipercaya, \[\begin{aligned} \hat{\beta}_0 &\sim N\left( \beta_0, \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}\right) \right) \\ \hat{\beta}_1 &\sim N\left( \beta_1, \frac{\sigma^2}{S_{xx}} \right). \end{aligned}\]

Dengan menstandarisasi hasil tersebut diperoleh \[\frac{\hat{\beta}_0 - \beta_0}{\text{SD}[\hat{\beta}_0]} \sim N(0, 1)\] dan \[\frac{\hat{\beta}_1 - \beta_1}{\text{SD}[\hat{\beta}_1]} \sim N(0, 1)\] dimana \[\text{SD}[\hat{\beta}_0] = \sigma\sqrt{\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}}\] dan \[\text{SD}[\hat{\beta}_1] = \frac{\sigma}{\sqrt{S_{xx}}}.\]

Karena tidak ada \(\sigma\) dalam latihan, maka kita akan mengestimasi menggunakan \(s_e\), dimasukkan pada persamaan yang dimiliki sebagai standar deviasi dari estimasi.

Kedua persamaan baru itu disebut kesalahan baku yang berupa estimasi standar deviasi dari distribusi sampling. \[\text{SE}[\hat{\beta}_0] = s_e\sqrt{\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}}\] \[\text{SE}[\hat{\beta}_1] = \frac{s_e}{\sqrt{S_{xx}}}\] Sekarang dibagi dengan kesalahan baku dan bukan standar deviasinya, hasil berikut digunakan untuk membuat selang kepercayaan dan melakukan uji hipotesis.

\[\frac{\hat{\beta}_0 - \beta_0}{\text{SE}[\hat{\beta}_0]} \sim t_{n-2}\] \[\frac{\hat{\beta}_1 - \beta_1}{\text{SE}[\hat{\beta}_1]} \sim t_{n-2}\] Untuk mendapatkan itu, perlu dicatat bahwa \[\frac{\text{RSS}}{\sigma^2} = \frac{(n-2)s_e^2}{\sigma^2} \sim \chi_{n-2}^2.\] dan ingat kembali variabel acak \(T\), \[T = \frac{Z}{\sqrt{\frac{\chi_{d}^2}{d}}}\] dilanjutkan distribusi \(t\) dengan derajat kebebasan \(d\), dimana \(\chi_d^2\) dan \(\chi^2\) variabel acak dengan derajat kebebasan \(d\).

Dituliskan, \[T\sim t_d\] dikatakan bahwa variabel acak \(T\) berhubungan dengan distribusi \(t\) dengan derajat kebebasan \(d\).

Lalu menggunakan trik klasik yaitu “dikalikan dengan 1” dan diatur ulang sehingga menghasilkan \[\begin{aligned} \frac{\hat{\beta}_1 - \beta_1}{\text{SE}[\hat{\beta}_1]} &= \frac{\hat{\beta}_1 - \beta_1}{s_e / \sqrt{S_{xx}}} \\ &= \frac{\hat{\beta}_1 - \beta_1}{s_e / \sqrt{S_{xx}}} \cdot \frac{\sigma / \sqrt{S_{xx}}}{\sigma / \sqrt{S_{xx}}} \\ &= \frac{\hat{\beta}_1 - \beta_1}{\sigma / \sqrt{S_{xx}}} \cdot \frac{\sigma / \sqrt{S_{xx}}}{s_e / \sqrt{S_{xx}}} \\ &= \frac{\hat{\beta}_1 - \beta_1}{\sigma / \sqrt{S_{xx}}} \bigg/ \sqrt{\frac{s_e^2}{\sigma^2}} \\ &= \frac{\hat{\beta}_1 - \beta_1}{\text{SD}[\hat{\beta}_1]} \bigg/ \sqrt{\frac{\frac{(n - 2)s_e^2}{\sigma^2}}{n - 2}} \sim \frac{Z}{\sqrt{\frac{\chi_{n-2}^2}{n-2}}} \sim t_{n-2} \end{aligned}\]

dimana \(Z\sim N(0,1)\).

Mengingat kembali distribusi \(t\) mirip dengan normal standar, tapi dengan ekor/ujung yang berbobot. Ketika derajat kebebasan meningkat, distribusi \(t\) semakin menjadi normal. Dibawah ini divisualisasikan distribusi normal standar seperti 2 contoh lainnya dari distribusi \(t\) dengan derajat kebebasan yang berbeda. Perhatikan bagaimana distribusi \(t\) dengan derajat kebebasan yang lebih besar mirip dengan kurva normal standar.

# define grid of x values
x = seq(-4, 4, length = 100)

# plot curve for standard normal
plot(x, dnorm(x), type = "l", lty = 1, lwd = 2,
     xlab = "x", ylab = "Density", main = "Normal vs t Distributions")
# add curves for t distributions
lines(x, dt(x, df = 1), lty = 3, lwd = 2, col = "darkorange")
lines(x, dt(x, df = 10), lty = 2, lwd = 2, col = "dodgerblue")

# add legend
legend("topright", title = "Distributions",
       legend = c("t, df = 1", "t, df = 10", "Standard Normal"), 
       lwd = 2, lty = c(3, 2, 1), col = c("darkorange", "dodgerblue", "black"))

1.4 Selang Kepercayaan untuk Kemiringan dan Perpotongan

Ingat kembali dimana selang kepercayaan berbentuk: \[\text{EST} \pm \text{CRIT} \cdot \text{SE}\] atau \[\text{EST} \pm \text{MARGIN}\] dimana \(\text{EST}\) sebuah estimasi untuk parameter ketertarikan. \(\text{SE}\) sebagai kesalahan baku dari suatu estimasi, dan \(\text{MARGIN}=\text{CRIT} \cdot \text{SE}\)

Lalu untuk \(\beta_0\) dan \(\beta_1\) dapat membuat selang kepercayaan menggunakan \[\hat{\beta}_0 \pm t_{\alpha/2, n - 2} \cdot \text{SE}[\hat{\beta}_0] \quad \quad \quad \hat{\beta}_0 \pm t_{\alpha/2, n - 2} \cdot s_e\sqrt{\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}}\] dan \[\hat{\beta}_1 \pm t_{\alpha/2, n - 2} \cdot \text{SE}[\hat{\beta}_1] \quad \quad \quad \hat{\beta}_1 \pm t_{\alpha/2, n - 2} \cdot \frac{s_e}{\sqrt{S_{xx}}}\] dimana \(t_{\alpha/2,n-2}\) adalah nilai kritis yang berarti \[P(t_{n-2} > t_{\alpha/2, n - 2}) = \alpha/2\].

1.5 Uji Hipotesis

Ingat kembali dimana uji statistik \(\text{(TS)}\) untuk uji rata-rata dalam bentuk: \[\text{TS} = \frac{\text{EST} - \text{HYP}}{\text{SE}}\] dimana \(\text{EsT}\) adalah estimasi untuk parameter ketertarikan, \(\text{HYP}\) adalah nilai hipotesis dari parameter, dan \(\text{SE}\) berupa kesalahan baku.

Maka ujinya \[H_0: \beta_0 = \beta_{00} \quad \text{vs} \quad H_1: \beta_0 \neq \beta_{00}\] gunakan uji statistiknya \[t = \frac{\hat{\beta}_0 - \beta_{00}}{\text{SE}[\hat{\beta}_0]} = \frac{\hat{\beta}_0-\beta_{00}}{s_e\sqrt{\frac{1}{n} + \frac{\bar{x}^2}{S_{xx}}}}\] yang berarti dibawah hipotesis 0, diikuti dengan distribusi \(t\) dengan \(n-2\) derajat kebebasan. Digunakan \(\beta_{00}\) untuk dinotasikan dengan nilai hipotesis \(\beta_0\).

Ujinya \[H_0: \beta_1 = \beta_{10} \quad \text{vs} \quad H_1: \beta_1 \neq \beta_{10}\] gunakan uji statistiknya \[t = \frac{\hat{\beta}_1-\beta_{10}}{\text{SE}[\hat{\beta}_1]} = \frac{\hat{\beta}_1-\beta_{10}}{s_e / \sqrt{S_{xx}}}s\] yang berarti dibawah hipotesis 0, diikuti dengan distribusi \(t\) dengan \(n-2\) derajat kebebasan. Digunakan \(\beta_{10}\) untuk dinotasikan dengan nilai hipotesis \(\beta_1\).

1.6 Contoh cars

Sekarang menggunakan cars dari bagian sebelumnya untuk mengilustrasi konsepnya. Awalnya gunakan model lm() lalu gunakan summary() untuk melihat hasilnya dengan lebih detail.

stop_dist_model = lm(dist ~ speed, data = cars)
summary(stop_dist_model)

## 
## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -29.069  -9.525  -2.272   9.215  43.201 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
## speed         3.9324     0.4155   9.464 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.38 on 48 degrees of freedom
## Multiple R-squared:  0.6511, Adjusted R-squared:  0.6438 
## F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

1.6.1 Uji dalam R

sekarang membahas hasil yang ditampilkan yang disebut Coefficient. Pertama ingat kembali bahwa dapat mengekstrak informasi secara langsung.

names(summary(stop_dist_model))

##  [1] "call"          "terms"         "residuals"     "coefficients" 
##  [5] "aliased"       "sigma"         "df"            "r.squared"    
##  [9] "adj.r.squared" "fstatistic"    "cov.unscaled"

summary(stop_dist_model)$coefficients

##               Estimate Std. Error   t value     Pr(>|t|)
## (Intercept) -17.579095  6.7584402 -2.601058 1.231882e-02
## speed         3.932409  0.4155128  9.463990 1.489836e-12

Fungsi names() memberi informasi yang diperoleh lalu gunakan operator $ dan Coefficients untuk mendapatkan informasi yang dibutuhkan. 2 hasil tersebut seharusnya terlihat tidak asing. \[\hat{\beta}_0 = -17.5790949\] dan \[\hat{\beta}_1 = 3.9324088\] yang berupa estimasi dari model parameter \(\beta_0\) dan \(\beta_1\).

Sekarang fokus pada baris kedua dari hasilm yang relevan terhadap \(\beta_1\).

summary(stop_dist_model)$coefficients[2,]

##     Estimate   Std. Error      t value     Pr(>|t|) 
## 3.932409e+00 4.155128e-01 9.463990e+00 1.489836e-12

Nilai pertama, Estimasi \[\hat{\beta}_1 = 3.9324088.\] Nilai kedua, Std. Error, kesalahan baku dari \(\hat{\beta_1}\), \[\text{SE}[\hat{\beta}_1] = \frac{s_e}{\sqrt{S_{xx}}} = 0.4155128.\] Nilai ketiga, t value, nilai uji statistik dari untuk menguji \(H_0: \beta_1 = 0\) vs \(H_1: \beta_1 \ne 0\),

\[t = \frac{\hat{\beta}_1-0}{\text{SE}[\hat{\beta}_1]} = \frac{\hat{\beta}_1-0}{s_e / \sqrt{S_{xx}}} = 9.46399.\] Terakhir, Pr(>|t|), memberikan nilai p-value dari uji tersebut. \[\text{p-value} = 1.4898365\times 10^{-12}\] Perlu dicatat disini secara spesifik untuk menguji apakah \(\beta_1=0\). Pada baris pertama dari laporan output dengan uji sama tetapi untuk \(\beta_0\).

summary(stop_dist_model)$coefficients[1,]

##     Estimate   Std. Error      t value     Pr(>|t|) 
## -17.57909489   6.75844017  -2.60105800   0.01231882

Kesimpulan, kode berikut menyimpan informasi summary(stop_dist_model)$coefficient dalam variabel baru stop_dist_model_test_info, lalu ekstrak tiap elemen pada variabel baru yang mendeskripsi informasi di dalamnya.

stop_dist_model_test_info = summary(stop_dist_model)$coefficients

beta_0_hat      = stop_dist_model_test_info[1, 1] # Estimate
beta_0_hat_se   = stop_dist_model_test_info[1, 2] # Std. Error
beta_0_hat_t    = stop_dist_model_test_info[1, 3] # t value
beta_0_hat_pval = stop_dist_model_test_info[1, 4] # Pr(>|t|)

beta_1_hat      = stop_dist_model_test_info[2, 1] # Estimate
beta_1_hat_se   = stop_dist_model_test_info[2, 2] # Std. Error
beta_1_hat_t    = stop_dist_model_test_info[2, 3] # t value
beta_1_hat_pval = stop_dist_model_test_info[2, 4] # Pr(>|t|)

Bisa dibuktikan beberapa persamaan yang setara: Uji \(t\) statistik \(\hat{\beta_1}\) dan 2 sisi nilai p / p-value berasosiasi dengan uji statistik.

(beta_1_hat - 0) / beta_1_hat_se

## [1] 9.46399

beta_1_hat_t

## [1] 9.46399

2 * pt(abs(beta_1_hat_t), df = length(resid(stop_dist_model)) - 2, lower.tail = FALSE)

## [1] 1.489836e-12

beta_1_hat_pval

## [1] 1.489836e-12

1.6.2 Signifikansi Regresi, Uji-t

Sekarang membahas uji signifikansi regresi. Pertama, perlu dicatat berdasarkan hasil distribusi, \(\beta_0\) dan \(\beta_1\) dapat diuji untuk nilai tertentu dan melakukan uji 1 dan 2 sisi.

Namun, satu uji spesifik, \[H_0: \beta_1 = 0 \quad \text{vs} \quad H_1: \beta_1 \neq 0\] sering digunakan. Mari pikirkan tentang pengujian ini dalam kaitannya dengan model regresi linier sederhana, \[Y_i = \beta_0 + \beta_1 x_i + \epsilon_i.\] Jika diasumsikan hipotesis 0 adalah benar, maka \(\beta_1=0\) dan memperoleh model, \[Y_i = \beta_0 + \epsilon_i.\] Dalam model ini, responden tidak bergantung dengan prediktor. Dengan begitu dapat dipikirkan dengan uji sebagai berikut,

• Dibawah \(H_0\) tidak ada perbedaan linier yang signifikan antara \(x\) dan \(y\).
• Dibawah \(H_1\) ada perbedaan linier yang signifikan antara \(x\) dan \(y\).

Untuk contoh cars,

• Dibawah \(H_0\) tidak ada perbedaan linier yang signifikan antara kecepatan dan jarak berhenti.
• Dibawah \(H_1\) ada perbedaan linier yang signifikan antara kecepatan dan jarak berhenti.

Uji tersebut dapat dilihat dari hasil summary(), \[\text{p-value} = 1.4898365\times 10^{-12}.\] Dengan nilai p yang sangat kecil, hipotesis 0 akan ditolak dalam batas wajar \(\alpha\), contoh \(\alpha=0.01\). Dengan begitu ada perbedaan linier yang signifikan antara kecepatan dan jarak berhenti.

Daritadi kami menekankan Linier.

# knitr::include_graphics("C:\Users\Acer\Desktop\Photoshop Project\Semester 5\Statistic Computation\Week 8\1.png")
library(OpenImageR)
img1 = readImage("1.png")
imageShow(img1)

Dalam plot dari data simulasi ini, dapat dilihat dengan jelas hubungan \(x\) dan \(y\), namun itu bukanlah hubungan linier. Jika kita menggambarkan garis pada data ini akan menjadi rata. Hasil uji untuk \(H_0:\beta_1=0\) vs \(H_1:\beta_1\ne0\) memberikan nilai p yang tinggi yaitu 0.7564548, maka gagal menolak dan dikatakan bahwa tidak ada hubungan linier antara \(x\) dan \(y\). Nantinya akan dilihat bagaimana memasukkan kurva pada data menggunakan model “linier”, untuk sekarang \(H_0:\beta_1=0\) vs \(H_1:\beta_1\ne0\) hanya dapat mendeteksi sebuah garis hubungan.

1.6.3 Selang Kepercayaan di R

Dengan menggunakan R, selang kepercayaan \(\beta_0\) dan \(\beta_1\) sangatlah mudah untuk didapatkan.

confint(stop_dist_model, level = 0.99)

##                  0.5 %    99.5 %
## (Intercept) -35.706610 0.5484205
## speed         2.817919 5.0468988

Dengan ini akan menghitung 99% selang kepercayaan dari \(\beta_0\) dan \(\beta_1\), baris pertama untuk \(\beta_0\) dan baris kedua untuk \(\beta_1\).

Untuk contoh cars saat menafsirkan dengan interval, dapat dikatakan 99% kepercayaan dalam penambahan kecepatan 1 mil per jam, rata-rata meningkat pada jarak berhenti diantara 2.8179187-5.0468988 kaki yang mana berupa selang untuk \(\beta_1\).

Perlu dicatat pada 99% selang kepercayaan tidak mengandung hipotesis 0. Karena tidak mengandung 0 yang artinya menolah uji \(H_0:\beta_1 = 0\) vs \(H_1:\beta_1\ne0\) saat \(\alpha=0.01\), seperti yang terjadi sebelumnya.

Anda harus agak curiga terhadap interval kepercayaan untuk \(\beta_0\), karena mencakup nilai negatif, yang sesuai dengan jarak berhenti negatif. Secara teknis interpretasinya adalah bahwa kami 99% yakin bahwa jarak berhenti rata-rata sebuah mobil yang menempuh jarak 0 mil per jam adalah antara -35.7066103 dan 0.5484205 kaki, tetapi kami tidak begitu percaya itu, karena kami benar-benar yakin bahwa itu akan terjadi. non-negatif.

Catatan, nilai spesifik dapat diperoleh dari output/hasil beberapa cara. Kode ini tidak run tetapi sebaiknya diperiksa hubungan output dari kode diatas.

confint(stop_dist_model, level = 0.99)[1,]

##       0.5 %      99.5 % 
## -35.7066103   0.5484205

confint(stop_dist_model, level = 0.99)[1, 1]

## [1] -35.70661

confint(stop_dist_model, level = 0.99)[1, 2]

## [1] 0.5484205

confint(stop_dist_model, parm = "(Intercept)", level = 0.99)

##                 0.5 %    99.5 %
## (Intercept) -35.70661 0.5484205

confint(stop_dist_model, level = 0.99)[2,]

##    0.5 %   99.5 % 
## 2.817919 5.046899

confint(stop_dist_model, level = 0.99)[2, 1]

## [1] 2.817919

confint(stop_dist_model, level = 0.99)[2, 2]

## [1] 5.046899

confint(stop_dist_model, parm = "speed", level = 0.99)

##          0.5 %   99.5 %
## speed 2.817919 5.046899

Bisa dipastikan juga kalkulasi selang \(\beta_1\) yang dilakukan dengan menggunakan R.

# store estimate
beta_1_hat = coef(stop_dist_model)[2]

# store standard error
beta_1_hat_se = summary(stop_dist_model)$coefficients[2, 2]

# calculate critical value for two-sided 99% CI
crit = qt(0.995, df = length(resid(stop_dist_model)) - 2)

# est - margin, est + margin
c(beta_1_hat - crit * beta_1_hat_se, beta_1_hat + crit * beta_1_hat_se)

##    speed    speed 
## 2.817919 5.046899

1.7 Interval Kepercayaan diri untuk Rata-Rata Respon

Selain interval kepercayaan untuk \(\beta_0\) dan \(\beta_1\), ada dua perkiraan interval umum lainnya yang digunakan dalam regresi. Yang pertama disebut interval kepercayaan untuk rata-rata respon . Seringkali, kami ingin perkiraan interval untuk rata-rata, \(E[Y|X=x]\) untuk nilai \(x\) tertentu. Dalam keadaan ini kita menggunakan \(\hat{y}(x)\) sebagai perkiraan dari \(E[Y|X=x]\). Kita mengubah sedikit notasi untuk memperjelas bahawa nilai prediksi adalah fungsi dari nilai \(x\)

\[\hat{y}=\hat\beta_0 + \hat\beta_{1}x\]

ingat, \[E[Y|X=x]=\beta_0 +\beta_{1}x\]

Jadi, \(\hat{y}\) adalah perkiraan yang bagus karena tidak bias \[E[\hat{y}]=\beta_0+\beta_{1}x\]

kemudiaan bisa mendapatkan, \[{Var[\hat{y}(x)]= σ^2({1\over n}+{x-\bar{x}^2 \over S_{xx}})}\]

Seperti perkiraan lainnya \(\hat{y}\) juga mengikuti distribusi normal. Sejak \(\hat\beta_0\) dan \(\hat\beta_1\) adalah kombinasi linear dari variabel random normal, \(\hat{y}\). \[\hat{y}∼N(β_0+β_1x,σ^2({1\over n}+{x-\bar{x}^2 \over S_{xx}}))\]

Dan akhirnya, sejak kita membutuhkan memperkirakan variansi, kita harus memperkirakan standar error. \[SE[\hat{y}]=s_e\sqrt{\frac{1}{n}+\frac{(x-x^2)}{S_{xx}}}\]

lalu kita dapat menggunakan untuk menemukan interval kepercayaan untuk rata-rata respon. \[\hat{y}±t_{α/2,n-2}. s_e\sqrt{\frac{1}{n}+\frac{(x-x^2)}{S_{xx}}}\]

Untuk menemukan interval kepercayaan untuk respons rata-rata menggunakan R, kami menggunakan fungsi predict(). Kami memberikan fungsi model yang dipasang serta data baru, disimpan sebagai bingkai data. (Ini penting, sehingga R mengetahui nama variabel prediktor.) Di sini, kami menemukan interval kepercayaan diri untuk jarak henti rata-rata ketika mobil bepergian 5 mil per jam dan ketika mobil bepergian 21 mil per jam.

new_speeds = data.frame(speed = c(5, 21))
predict(stop_dist_model, newdata = new_speeds, 
        interval = c("confidence"), level = 0.99)

##         fit       lwr      upr
## 1  2.082949 -10.89309 15.05898
## 2 65.001489  56.45836 73.54462

1.8 Prediksi Interval untuk Observasi Baru

Terkadang kita ingin meng-estimasikan interval yang baru, \(Y\), untuk nilai tertentu dari \(x\). Ini sangat mirip dengan interval respon rata-rata, \(E[Y | X = x]\), tapi berbeda di satu hal yang sangat penting.

Tebakan terbaik kita untuk observasi yang baru masih tetap \(\hat y(x)\). Estimasi rata-rata masih tetap prediksi terbaik yang bisa kita buat. Perbedaanya terdapat pada jumlah variabilitasnya. Kita tahu bahwa observasi itu akan bervariasi tentang garis regresi yang sebenarnya menurut distribusi \(N(0,σ^2)\). Karena ini kita menambahkan beberapa faktor dari \(σ^2\) untuk estimasi variabilitas kita dengan tujuan untuk menjelaskan variabilitas dari observasi tentang garis regresi.

\[Var[\hat y(x)+ϵ]= Var[\hat y(x)]+Var [ϵ]\] \[= \sigma^2 (\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}}) + σ^2 \] \[ = \sigma^2 (1 +\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}})\] \[\hat y(x) + ϵ ∼ N (\beta_0 + \beta_1x,σ^2(1 +\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}}))\] \[SE[\hat y(x)+ϵ] = s_e \sqrt {1 +\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}}}\]

Lalu kita dapat menemukan interval prediksi menggunakan,

\[\hat y(x) ± t_{\alpha/2,n-2} ⋅ s_e\sqrt{1 +\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}}}.\] Untuk menghitung ini untuk kumpulan poin di R, ketahuilah bahwa hanya ada sedikit perubahan didalam sintaks dari penemuan interval kepercayaan untuk respon rata-rata.

predict(stop_dist_model, newdata = new_speeds, 
        interval = c("prediction"), level = 0.99)

##         fit       lwr       upr
## 1  2.082949 -41.16099  45.32689
## 2 65.001489  22.87494 107.12803

Juga ketahuilah bahwa dua interval ini lebih lebar dibanding interval kepercayaan yang sesuai untuk respon rata-rata.

1.9 Keyakinan dan Band Prediksi

Seringkali kita akan memplot kedua interval kepercayaan untuk respon rata-rata dan interval prediksi untuk semua kemungkinan dari nilai \(x\). Kita memanggil ini sebagai kepercayaan dan band prediksi.

speed_grid = seq(min(cars$speed), max(cars$speed), by = 0.01)
dist_ci_band = predict(stop_dist_model, 
                       newdata = data.frame(speed = speed_grid), 
                       interval = "confidence", level = 0.99)
dist_pi_band = predict(stop_dist_model, 
                       newdata = data.frame(speed = speed_grid), 
                       interval = "prediction", level = 0.99) 

plot(dist ~ speed, data = cars,
     xlab = "Speed (in Miles Per Hour)",
     ylab = "Stopping Distance (in Feet)",
     main = "Stopping Distance vs Speed",
     pch  = 20,
     cex  = 2,
     col  = "grey",
     ylim = c(min(dist_pi_band), max(dist_pi_band)))
abline(stop_dist_model, lwd = 5, col = "darkorange")

lines(speed_grid, dist_ci_band[,"lwr"], col = "dodgerblue", lwd = 3, lty = 2)
lines(speed_grid, dist_ci_band[,"upr"], col = "dodgerblue", lwd = 3, lty = 2)
lines(speed_grid, dist_pi_band[,"lwr"], col = "dodgerblue", lwd = 3, lty = 3)
lines(speed_grid, dist_pi_band[,"upr"], col = "dodgerblue", lwd = 3, lty = 3)
points(mean(cars$speed), mean(cars$dist), pch = "+", cex = 3)

Beberapa hal yang perlu diperhatikan: ~Kita menggunakan argumen ylim untuk melonggarkan sumbu y dari plot tersebut, sejak bandnya lebih jauh daripada titik yang ada. ~Kita menambahkan poin di poin \((x,\bar y)\) Ini merupakan sebuah titik dimana garis regresinya akan selalu melewatinya. (Pikirkan kenapa bisa terjadi.) Ini adalah titik dimana kedua kepercayaan dan band prediksi menjadi yang tersempit. Lihatlah kepada standar kesalahannya dari keduanya untuk mengetahui kenapa itu bisa terjadi.

1.10 Uji F, Signifikansi Regresi

Kasus dari regresi linear sederhana, uji \(t\) untuk signifikasi regresi merupakan nilai equivalent untuk uji lainnya, yaitu uji \(F\) untuk signifikansi regresi. Equivalent ini hanya akan menjadi benar untuk regresi linear sederhana, dan di bagian selanjutnya kita hanya akan menggunakan uji \(F\) untuk signifikansi regresinya.

Ingatlah kembali di seksi terakhir dari dekomposisi varians yang kita lihat sebelum menghitung nilai `\(R^2\),

\[∑_{i=1}^n ( y_i - \bar y)^2 = \Sigma_{i-1}^n (y_i - \bar y_i)^2 + \Sigma_{i=1}^n (\hat y_i - \bar y)^2,\]

atau, lebih singkatnya,

\[SST = SSE + SSReg.\]

Untuk mengembangkan uji \(F\), kita akan mengatur lagi informasi ini di tabel ANOVA,

# knitr::include_graphics("C:\Users\Acer\Desktop\Photoshop Project\Semester 5\Statistic Computation\Week 8\2.png")
library(OpenImageR)
img2 = readImage("2.png")
imageShow(img2)

ANOVA atau analisis dari varians akan menjadi sebuah konsep yang sering kita bahas dalam kelas ini. Untuk sekarang, kita akan fokus kepada hasil-hasil didalam tabel, dimana \(F\) statistiknya,

\[\frac {\Sigma_{i=1}^n ( y_i - \bar y)^2/1}{\Sigma_{i=1}^n ( y_i - \hat y_i)^2/(n-2)} ∼ F_{1,n-2}\]

dimana mengikuti distribusi \(F\) dengan derajat kebebasan 1 dan \(n-2\) dibawah hipotesis nol. Sebuah distribusi \(F\) adalah distribusi kontinu yang hanya mengambil nilai positif dan mempunyai dua parameter, yang merupakan dua derajat kebebasan.

Ingatlah kembali, didalam signifikansi dari uji regresi, \(Y\) tidak bergantung kepada \(x\) di hipotesis nol.

\[H_0 : \beta_1 = 0\] \[ Y_i = \beta_0 + \epsilon_i\]

Disaat alternatif dari hipotesis \(Y\) bergantung kepada \(x\).

\[H_0 : \beta_1 ≠ 0\] \[ Y_i = \beta_1x_i + \epsilon_i\]

Kita dapat menggunakan \(F\) statistik untuk menunjukkan pengujian ini.

\[F = \frac {\Sigma_{i=1}^n (\hat y_i - \bar y)^2/1}{\Sigma_{i=1}^n ( y_i - \hat y_i)^2/(n-2)}\]

Secara khusus, kita akan menolak hipoteis nolnya disaat nilai \(F\) statistiknya besar, itu dia, disaat ada probabilitas rendah bahwa pengamatan secara kebetulan, dapat berasal dari model hipotesis nol. Kita akan membiarkan R menghitung nilai p untuk kita.

UNtuk menunjukkan uji \(F\) di R, anda dapat melihat di baris terakhir dari data yang keluar dari summary() yang dipanggil F-statistic yang memberikan nilai dari uji statistik, nilai derajat kebebasan yang relevan, juga nilai p dari pengujian tersebut.

summary(stop_dist_model)

## 
## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -29.069  -9.525  -2.272   9.215  43.201 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
## speed         3.9324     0.4155   9.464 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.38 on 48 degrees of freedom
## Multiple R-squared:  0.6511, Adjusted R-squared:  0.6438 
## F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

Selain itu, anda dapat menggunakan fungsi anova() untuk menampilkan informasi dari tabel ANOVA.

anova(stop_dist_model)

## Analysis of Variance Table
## 
## Response: dist
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## speed      1  21186 21185.5  89.567 1.49e-12 ***
## Residuals 48  11354   236.5                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Ini juga memberikan nilai p dari suatu pengujian. Anda harus memperhatikan bahwa nilai p dari uji \(t\) itu sama. Anda juga mungkin sadar bahwa nilai dari uji statistik dari uji \(t\), 9.46399, bisa di kuadratkan untuk mendapatkan nilai dari \(F\) statistik, 89.5671065.

Perhatikan bahwa ada beberapa jalanlain yang sama yang dapat dilakukan di R, yang sering kita lohat untuk membandingkan dua model.

anova(lm(dist ~ 1, data = cars), lm(dist ~ speed, data = cars))

## Analysis of Variance Table
## 
## Model 1: dist ~ 1
## Model 2: dist ~ speed
##   Res.Df   RSS Df Sum of Sq      F   Pr(>F)    
## 1     49 32539                                 
## 2     48 11354  1     21186 89.567 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pernyataan model lm(dist ~ 1, data = cars) menerapkan model \(Y_i = \beta_0 + \epsilon_i\) kepada data cars. Perhatikan bahwa \(\hat y = \bar y\) disaat \(Y_i = \beta_0 + \epsilon_i\).

Pernyataan model lm(dist ~ speed, data = cars) menerapkan model \(Y_i = \beta_0 + \beta_1x_i + \epsilon_i\).

Lalu kita dapat memikirkan kegunaan dari anova() sebagai perbandingan secara langsung dari dua model. (Perhatikan bahwa kita mendapatkan lagi nilai p yang sama.)