CS-Lab8-Inference-in-Linear-Regression

1 Inferensial di Regresi Linear

Setelah membaca seksi ini, kamu akan bisa:

~Mengerti mengenai distribusi dari suatu estimasi regresi.

~Membuat interval estimasi untuk regresi parameter, respon rata-rata, dan prediksi.

~Menguji signifikan regresinya.

Di seksi terakhir, kita memecahkan model regresi linear sederhana, yaitu

\[Y_i = β_0 + β_1x_i + ϵ_i\] dimana \(ϵ_i ∼ N(0,σ^2)\). Lalu kita menggunakan observasi \((x_i,y_i)\), untuk \(i = 1,2,...,n\), untuk menemukan ilai dari \(β_0 dan β_1\) yang bisa diminimumkan

\[f(β_0,β_1)=∑_{i=1}^n(y_i-(β_0 + β_1x_i))^2\].

Kita memanggil nilai \(\hatβ_0\) dan \(\hatβ_1\), dimana kita menemukannya menjadi

\[\hatβ_1= \frac {S_{xy}}{S_{xx}} = \frac {∑_{i=1}^n(x_i - \bar x)(y_i - \bar y )}{∑_{i=1}^n(x_i-x)^2}\] \[\hatβ_0 = \bar y - \hatβ_1 \bar x\].

Kita juga akan memperkirakan nilai \(σ^2\) menggunakan \(s_e^2\). Di kata lain, kita menemukan bahwa \(s_e\) adalah sebuah estimasi untuk \(σ\), dimana

\[s_e = RSE = \sqrt \frac {1}{n - 2} ∑_{i=1}^n e_i^2\] dimana kita juga memanggil RSE, sebagai “Residual Standart Error.”

Disaat mengaplikasikan ke data cars, kita mendapatkan hasil seperti berikut:

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

Di seksi terakhir, kita hanya mendiskusikan tentang nilai estimate, Residual standard error, dan Multiple R-squared. Di seksi ini, kita akan mendiskusikan semua informasi dibawah nilai Coefficients, begitu juga dengan nilai F-statistic.

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")

Dari grafik diatas, kita dapat menyimpulkan bahwa ada suatu korelasi yang ada diantara jarak berhenti dengan suatu kecepatan, dimana hasilnya adalah variabel kecepatan, mempengaruhi jarak berhenti, semakin cepat kecepatan pada sesuatu, maka semakin jauh pula jarak berhentinya.

Untuk memulai, kita akan perhatikan bahwa ada ekspresi equivalen yang lain untuk \(S_{xy}\) dimana tidak kita lihat di sesi yang terakhir, \[S_{xy} = ∑_{i=1}^n (x_i - \bar x)(y_i - \bar y) = ∑_{i=1}^n (x_i - \bar x)y_i\].

Ini mungkin akan menjadi nilai equivalen yang mengejutkan. (Mungkin akan mencoba untuk membuktikannya.) Lagipula, itu akan menjadi berguna sebagai sebuah konsep ilustrasi di seksi ini.

Perhatikan bahwa, \(\hat β_1\) adalah sampel statistik yang sudah terhitung dengan suatu observasi data seperti yang dituliskan diatas, yaitu \(\hat β_0\).

Bagaimanapun juga, seksi ini akan sering memudahkan kita untuk menggunakan kedua \(\hat β_1\) dan \(\hat β_0\) sebagai variabel acak, itu dia, kita masih belum mempunyai nilai observasi untuk masing-masing dari\(Y_i\). Disaat kasusnya adalah ini, kita akan menggunakan notasi yang berbeda, mensubtitusikan \(Y_i\) untuk kasus \(y_i\) yang lebih rendah.

\[\hatβ_1= \frac {∑_{i=1}^n(x_i - \bar x)Y_i}{∑_{i=1}^n(x_i-x)^2}\] \[\hatβ_0 = \bar Y - \hatβ_1 \bar x\]

Di seksi terakhir kita menyetujui bahwa perkiraan ini dari model parameter \(β_0\) dan \(β_1\) yang tidak diketahui, itu bagus karena kita mendapatkan mereka dari meminimumkan kesalahan. Kita sekarang akan mendiskusikan tentang teorema Gauss-Markov yang mengambil lebih lanjut tentang ide ini, menunjukan bahwa perkiraan ini sebenarnya merupakan perkiraan yang “terbaik”, dari sudut pandang tertentu.

1.1 Teorema Gauss-Markov

Teorema Gauss-Markov ini memberitahukan kepada kita bahwa disaat kita memperkirakan parameter dari regresi linear sederhana model \(β_0\) dan \(β_1\), nilai \(β_0\) dan \(β_1\) yang kita percayakan adalah perkiraan linear tidak bias yang teraik, atau lebih singkatnya “BLUE / Best Linear Unbiased Estimates”. (Sebuah kondisi sebenarnya untuk teorema Gauss-Markov yang lebih rileks dibanding model SLR.)

Sekarang kita akan mendiskusikan tentang linear, tidak bias, dan terbaik, yang berhubungan dengan perkiraan tersebut.

1.1.1 Linear

Mengulang kembali, di pengaturan SLR dimana nilai \(x_i\) sudah diketahui, dan dianggap sebagai suatu kuantitas. Lalu sebuah perkiraan linear adalah salah satu yang dapat ditulis sebagai kombinasi linear dari \(Y_i\). Di sebuah kasus dari \(\hat β_1\) kita melihat

\[\hatβ_1= \frac {∑_{i=1}^n(x_i - \bar x)Y_i}{∑_{i=1}^n(x_i-x)^2} = ∑_{i=1}^n k_iY_i = k_1Y_1 + k_2Y_2 + ... + k_nY_n \]

dimana \(k_i = \frac {(x_i - \bar x)}{∑_{i=1}^n (x_i - \bar x)^2}\).

Di pertunjukan yang mirip, kita dapat menunjukkan bahwan \(\hatβ_0\) dapat ditulis sebagai kombinasi linear dari \(Y_i\). Jadi kedua \(\hatβ_0\) dan \(\hatβ_1\) adalah suatu estimator dari fungsi linear.

1.1.2 Tidak Bias

Sekarang kita sudah tau bahwa perkiraan kita adalah fungsi linear, seberapa baguskah estimasi ini? Salah satu pengukuran “kebaikan” dari suatu perkiraan adalah bias. Lebih spesifik, kita cenderung memperkirakan yang tidak bias, berarti ekspektasi nilai mereka adalah parameter yang sedang di perkirakan.

Di dalam kasus dari perkiraan regresi, kita memiliki,

\[E[\hatβ_0]=β_0\] \[E[\hatβ_1]=β_1.\]

Ini memberitahukan kepada kita bahwa, disaat kondisi dari suatu model SLR bertemu, rata-rata dari perkiraan kita akan benar. Lagipula, yang sudah kita lihat dari seksi terakhir dimana kita mensimulasikan dari model SLR, itu tidak berarti bahwa setiap perkiraan akan benar. Hanya itu, jika kita mengulang prosesnya ke jumlah yang tidak terhingga, perkiraan rata-rata itu akan benar.

1.1.3 Terbaik

Sekarang, jika kita membatasi diri kita kepada kedua fungsi linear dan perkiraan tidak bias, bagaimana kita akan memastikan perkiraan terbaiknya? Perkiraan dengan varians terminim.

Catatan pertama adalah ini sangat mudah untuk membuat suatu perkiraan untuk \(β_1\) yang memiliki varians yang rendah, tapi bukan tidak bias. Sebagai contoh, pastikan:

\[\hatθ_{BAD}=5\]

Lalu, sejak \(\hatθ_{BAD}\) adalah sebuah nilai konstan,

\[Var[\hatθ_{BAD}]=0\]

Bagaimanapun juga sejak,

\[E[\hatθ_{BAD}] = 5\]

kita mengatakan bahwa \(\hatθ_{BAD}\) adalah estimator yang bias kecuali kalau \(β_1=5\), dimana sebelumnya kita tidak akan tau. Untuk alasan ini, ini merupakan perkiraan yang mengerikan (kecuali kalau ada kesempatan \(β_1=5\)) bahkan jika itu adalah varians terkecil yang memungkinkan. Ini adalah bagian dari suatu alasan yang membatasi diri kita untuk memperkirakan yang tidak bias. Apa yang bagus pada sebuah estimasi, jika itu memperkirakan kuantitas yang salah?

Jadi sekarang, untuk pertanyaan alaminya, apakah varians dari \(\hatβ_0\) dan \(\hatβ_1\)? Mereka adalah,

\[Var[\hatβ_0] = σ^2(\frac{1}{n} + \frac {\bar x^2}{S_{xx}})\] \[Var[\hatβ_1] = \frac {σ^2}{S_{xx}}\]

Ini mengukur variabilitas dari dua perkiraan ke kesempatan acak disaat sedang menyampling data. Apakah ini “yang terbaik”? Apakah varians ini adalah varians yang terkecil yang bisa kita dapatkan? Kamu hanya bisa percaya kepada kata-kata kami bahwa ini memang sudah benar.

1.2 Distribusi Sampling

Sekarang kita harus mendefinisikan ulang perkiraan untuk \(\hatβ_0\) dan \(\hatβ_1\) yang menjadi variabel acaknya, kita dapat berdiskusi tentang distribusi sampling mereka, yang mana distribusi ketika dianggap sebagai variabel acak.

Sejak kedua \(\hatβ_0\) dan \(\hatβ_1\) adalah kombinasi linear dari \(Y_i\) dan setiap \(Y_i\) adalah distribusi normal, maka kedua \(\hatβ_0\) dan \(\hatβ_1\) juga merupakan / mengikuti distribusi normal.

Lalu, menaruh semuanya yang diatas secara bersama, kita sampai di distribusi dari \(\hatβ_0\) dan \(\hatβ_1\).

Untuk \(\hatβ_1\), kita katakan

\[\hatβ_1 = \frac {S{xy}}{S{xx}} = \frac {∑_{i=1}^n(x_i - \bar x)Y_i}{∑_{i=1}^n(x_i-x)^2} ∼ N (β_1, \frac {σ^2}{∑_{i=1}^n(x_i-\bar x)^2})\]

Atau lebih singkatnya,

\[\hatβ_1 ∼ N (β_1, \frac {σ^2}{S{xx}})\]

Dan untuk \(\hatβ_0\),

\[\hatβ_0 = \bar Y - \hatβ_1 \bar x ∼ (β_0, \frac {σ^2 ∑_{i=1}^n x_i^2}{n∑_{i=1}^n (x_i-\bar x)^2})\]

Atau lebih singkatnya,

\[\hatβ_0 ∼ N (β_0, σ^2(\frac {1}{n}+\frac {\bar x^2}{S{xx}}))\]

Di bagian ini kita telah lalai untuk membuktikan nomor dari hasil tersebut. Alih-alih untuk mengerjakan dengan derivasi yang membosankan dari distribusi sampling ini, kita akan membenarkan hasil ini untuk diri kita sendiri dengan menggunakan simulasi.

Catatan untuk para pembaca: Derivasi dan pembuktian ini mungkin akan ditambahkan ke lampiran di suatu saat nanti. Kamu juga bisa menemukan hasil ini di waktu yg akan datang di buku teks tentang regresi linear sederhana. Di UIUC, hasil ini akan dipresentasikan di kedua STAT 424 dan STAT 425. Namun, sejak anda tidak akan bertanya untuk menampilkan derivasi dari tipe yang ada di kelas ini, mereka sekarang dihilangkan.

1.2.1 Simulasi Distribusi Sampling

Untuk memastikan hasil diatas, kita akan mensimulasikan sample dengan ukuran \(n = 100\) dari model

\[Y_i = β_0 + β_1 x_i + ϵ_i\]

dimana \(ϵ_i ∼ N(0,σ^2)\). Di kasus ini, parameternya dapat diketahui sebagai:

~ $β_0 = 3$
~ $β_1 = 6$
~ $σ^2 = 4$

Lalu, berdasarkan pernyataan diatas, kita harus menemukan bahwa

\[\hat β_1 ∼ N(β_1, \frac {σ^2}{S{xx}})\]

dan

\[\hat β_0 ∼ N(β_0,σ^2(\frac {1}{n}+\frac {\bar x^2}{S{xx}}))\]

Pertama kita sebelumnya harus menentukan berapa nilai x kita di simulasi ini, sejak nilai x di SLR juga dianggap sebagai kuantitas yang diketahui. Pilihan untuk nilai x adalah semaunya kita. Disini kita juga menentukan benih untuk pengacakan dan perhitungan \(S{xx}\) dimana kita akan tetap bergerak kedapan.

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

Kita juga memastikan nilai parameter kita.

beta_0 = 3
beta_1 = 6
sigma  = 2

Dengan informasi ini, kita tau distribusi samplingnya harus menjadi:

(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 β_1 ∼ N(6,0.1176238)\] dan \[\hat β_0 ∼ N(3,0.04).\]

Ini dia,

\[E[\hat β_1] = 6\] \[Var[\hat β_1] = 0.1176238\] dan \[E[\hat β_0] = 3\] \[Var[\hat β_0] = 0.04.\]

Kita sekarang mensimulasikan data dari model ini sebanyak 10,000 kali. Perhatikan bahwa ini mungkin bukanlah jalan R yang biasanya dilakukan disaat melakukan simulasi. Kita akan menunjukkan simulasi dengan cara ini yang berupaya untuk menjelaskan. Sebagai contoh, kita dapat menggunakan fungsi sim_slr() dari seksi sebelumnya. Kita juga akan menyimpan variabel ini di lingkungan global untuk membuat bingkai data untuk setiap simulasi data yang baru.

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]
}

Setiap kali kita mensimulasikan data, kita mendapatkan hasil dari koefisien perkiraannya. Variabel dari beta_0_hats dan beta_1_hats akan menyimpan 10,000 nilai simulasi masing-masing dari \(\hat β_0\) dan \(\hat β_1\).

Pertama kita akan memastikan distribusi dari \(\hat β_1\).

mean(beta_1_hats) # empirical mean

## [1] 6.001998

beta_1            # true mean

## [1] 6

var_beta_1_hat    # true variance

## [1] 0.1176238

Kita melihat bahwa cara dan rata-rata yang benar, dan varians yang sangat mirip. Kita juga memastikan bahwa distribusi empiris bersifat normal. Untuk melakukannya, kita membuat plot histogram dari beta_1_hats, dan menambahkan kurva untuk distribusi yang benar dari \(\hat β_1\). Kita menggunakan prob = TRUE untuk menaruh histogramnya di skala yang sama dengan kurva normal.

# 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 kita mengulang proses dari \(\hat β_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)

Di simulasi belajar ini, kita hanya mensimulasikan jumlah sampel yang terbatas. Untuk benar-benar memastikan hasil distribusinya, kita harus mengamati jumlah sampel yanhg tidak terbatas. Namun, plot ini harus menjelaskan bahwa jika kita melanjutkan melakukan simulasinya, hasil empirisnya akan mendekati apa yang kita harapkan.

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 Standar Kesalahan

Jadi sekarang kita percaya bahwa hasil dari dua distribusi,

\[\hat β_0 ∼ N (β_0,σ^2(\frac {1}{n}+\frac {\bar x^2}{S{xx}}))\] \[\hat β_1 ∼ N(β_1,\frac {σ^2}{S{xx}}).\]

Lalu dengan menyederhanakan hasil ini, kita menemukan bahwa

\[\frac {\hat β_0 - β_0}{SD[\hat β_0]} ∼ N(0,1)\]

dan \[\frac {\hat β_1 - β_1}{SD[\hat β_1]} ∼ N(0,1)\] dimana \[SD[\hat β_0] = σ \sqrt {\frac {1}{n} + \frac {\bar x^2}{S{xx}}}\] dan \[SD[\hat β_1] = \frac {σ}{\sqrt {S{xx}}}.\]

Sejak kita tidak tau nilai σ di latihan kita, kita akan memperkirakan itu menggunakan \(s_e\), yang mana kita hubungkan ke ekspresi kita yang ada untuk standar deviasi dari estimasi kita.

Kedua ekspresi baru ini dipanggil standar errors atau kesalahan yang merupakan perkiraan dari standar deviasi untuk distribusi samplingnya.

\[SE[\hat β_0] = s_e \sqrt {\frac {1}{n}+\frac {\bar x^2}{S{xx}}}\] \[SE[\hat β_1] = \frac {s_e}{\sqrt{S{xx}}}\]

Sekarang jika kita membagikan dengan standar errornya, alih-alih standar deviasinya, kita mendapatkan hasil dimana akan memungkinkan kita untuk membuat interval kepercayaan dan menunjukkan pengujian hipotesis.

\[\frac {\hat β_0 - β_0}{SE[\hat β_0]} ∼ t_{n-2}\]

Untuk melihat ini, pertama kita perhatikan bahwa,

\[\frac {RSS}{σ^2} = \frac {(n-2)s_e^2}{σ^2} ∼ χ_{n-2}^2\]

Juga memanggil ulang variabel acak \(T\) dengan,

\[T = \frac {Z}{\sqrt {\frac {χ_d^2}{d}}}\]

mengikuti distribusi \(t\) dengan derajat kebebasan \(d\), dimana \(χ_d^2\) adalah \(χ^2\) variabel acak dari derajat kebebasan \(d\).

Kita menulis, \[T ∼ t_d\]

untuk mengatakan bahwa variabel acak \(T\) mengikuti distribusi \(t\) dengan derajat kebebasan \(d\).

Lalu kita menggunakan trik klasik untuk “mengkalikan dengan 1” dan beberapa aturan ulang untuk sampai ke

\[ \frac {\hat β_1 -β_1}{SE[{\hat β_1]}} = \frac {\hat β_1 -β_1}{s_e/ \sqrt{S_{xx}}}\] \[= \frac {\hat β_1 -β_1}{s_e/ \sqrt{S_{xx}}} . \frac {σ/\sqrt{S_{xx}}}{σ/\sqrt{S_{xx}}}\] \[=\frac {\hat β_1 -β_1}{SD[\hat β_1]} / \sqrt {\frac {\frac {(n-2)s_e^2}{σ^2}}{n-2}} ∼ \frac {Z}{\sqrt \frac{χ_{n-2}^2}{n-2}} ∼ t_{n-2}\]

dimana \(Z∼N(0,1)\)

Ingatlah kembali bahwa distribusi \(t\) mirip dengan standar normal, tapi memiliki ekor yang lebih berat. Dengan derajat kebebasan yang semakin meningkat, distribusi \(t\) menjadi lebih dan lebih seperti standar normal. Dibawah ini, kita memplotkan distribusi standar normal sebagai dua contoh dari distribusi \(t\) dengan derajat kebebasan yang berbeda. Ketahuilah bagaimana distribusi \(t\) dengan derajat kebebasan yang lebih besar tersebut lebih mirip dengan kurva standar normal.

# 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 Interval Keyakinan untuk Slope and Intercept

Mengulang kembali bahwa interval keyakinan untuk rata-rata sering diambil dalam bentuk: \[EST ± CRIT⋅SE\] atau \[EST ± MARGIN\]

dimana EST adalah perkiraan untuk parameter yang diinginkan, SE adalah standar kesalahan untuk suatu perkiraan, dan \(MARGIN = CRIT⋅SE\).

Lalu, untuk \(β_0\) dan \(β_1\) kita dapat membuat interval kepercayanya dengan menggunakan

\[\hatβ_0 ± t_{α/2,n-2} ⋅ SE[\hatβ_0]\]
\[\hatβ_0 ± t_{α/2,n-2} ⋅s_e\sqrt{\frac {1}{n} + \frac {x^2}{S{xx}}}\] dan \[\hatβ_1 ± t_{α/2,n-2} ⋅ SE[\hatβ_1]\] \[\hatβ_0 ± t_{α/2,n-2} ⋅\frac {s_e}{\sqrt {S_{xx}}}\]

dimana $ t_{α/2,n-2}$ adalah nilai kritis demikian juga \(P(t_{n-2} > t_{α/2,n-2}) = α/2\).

1.5 Uji Hipotesis

“Kita bisa berbicara tentang hipotesis ini sebagai ‘null hypothesis’ atau hipotesis nol, dan perlu dicatat bahwa hipotesis nol tidak pernah dibuktikan atau ditetapkan, tetapi mungkin disangkal, selama eksperimen.” - Ronald Aylmer Fisher

Menginat kembali bahwa uji statistik \((TS)\) untuk menguji rata-rata sering diambil dengan bentuk:

\[TS = \frac {EST - HYP}{SE}\]

dimana EST adalah perkiraan untuk parameter yang diinginkan, HYP adalah nilai hipotesis untuk parameternya, dan SE adalah standar kesalahan untuk perkiraannya.

Jadi, untuk pengujian

\[H_0 : β_0 = β_{00} vs H_1 : β_0 ≠ β_{00}\]

kita menggunakan uji statistik

\[t = \frac {\hat β_0 - β_{00}}{SE[\hat β_0]} = \frac {\hat β_0 - β_{00}}{s_e \sqrt{\frac {1}{n} + \frac{\bar x^2}{S{xx}}}}\]

dimana, dibawah hipotesis nol, mengikuti distribusi \(t\) dengan derajat kebebasan \(n-2\). Kita menggunakan \(β_{00}\) untuk mencari tau nilai hipotesis untuk \(β_0\).

Begitupula, untuk pengujian

\[H_0 : β_1 = β_{10}\] vs \[H_1 : β_1 ≠ β_{10}\]

kita menggunakan uji statistik

\[t = \frac {\hat β_1 - β_{10}}{SE[\hat β_1]} = \frac {\hat β_1 - β_{10}}{s_e / \sqrt {S{xx}}}\]

dimana lagi, dibawah hipotesis nol, mengikuti distribusi \(t\) dengan derajat kebebasan \(n-2\). Kita sekarang menggunakan \(β_{10}\) untuk menentukan nilai hipotesis untuk \(β_1\)

1.6 Contoh cars

Kita sekarang kembali ke contoh cars dari seksi sebelumnya untuk mengilustrasikan konsep ini. Pertama kita mengepaskan dengan menggunakan model lm() lalu menggunakan summary() untuk melihat hasil dengan rincian yang lebih hebat.

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 Tes di R

Kita sekarang akan mendiskusikan tentang hasil yang ditampilkan yang dipanggil Coefficients. Pertama, mengulang lagi bahwa kita dapat meng-ekstrak informasi ini 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 dari names() memberitahukan kepada kita informasi apa yang sedang ada, dan lalu kita menggunaan operator $ dan coefficients untuk meng-ekstrak informasi yang sudah kita ingikan. Dua nilai ini harus saling dikenal.

\[\hatβ_0 = -17.5790949\] dan \[\hatβ_1 = 3.9324088\]

dimana merupakan perkiraan kita untuk model parameter \(β_0\) dan \(β_1\). Sekarang ayo fokus di baris kedua dari data yang keluar, dimana berhubungan dengan \(β_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

Lagi, jumlah nilai pertama, Estimate adalah

\[ \hat β_1 = 3.9324088\]

Nilai kedua, Std.Error, adalah nilai standar kesalahan dari \(\hatβ_1\),

\[SE[\hatβ_1] = \frac {s_e}{\sqrt {S{xx}}} = 0.4155128.\]

Nilai ketiga, t value, adalah nilai dariuji statistik untuk pengujian \(H_0 : β_1 = 0 vs H_1 : β_1 ≠ 0\),

\[t = \frac {\hat β_1 - 0}{SE[\hat β_1]} = \frac {\hat β_1 - 0}{s_e / \sqrt {S{xx}}} = 9.46399.\]

Terakhir, Pr(>|t|), memberikan kita nilai p dari pengujian tersebut.

\[p-value = 1.4898365 × 10^{-12}\]

Catat disini, kita secara khusus menguji iya atau tidaknya \(β_1 = 0\).

Baris pertama dari laporan data keluar nilainya sama, tetapi untuk \(β_0\).

summary(stop_dist_model)$coefficients[1,]

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

Dalam penjelasan, kode berikus menyimpan informasi dari summary(stop_dist_model)$coefficients di variabel baru stop_dist_model_test_info, lalu ekstrak setiap elemen menjadi variabel baru yang mendeskripsikan informasi yang dikandungnya.

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|)

Kita dapat memastikan beberapa ekspresi ekuivalen: uji \(t\) statistik untuk \(\hat β_1\) dan kedua sisi dari nilai p terkait dengan uji statistik tersebut.

(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 Uji-t, Signifikansi Regresi

Kita berhenti untuk berdiskusi tentang pentingnya uji regresi. Pertama, perhatikan bahwa berdasarkan dari hasil distribusi diatas, kita dapat menguji \(β_0\) dan \(β_1\) terhadap nilai apapun, dan menunjukkan pengujian untuk satu dan dua sisi.

Namun, satu uji yang paling spesifik,

\[H_0 :β_1 = 0\] vs \[H_1 : β_1 ≠ 0\] yang sangat sering digunakan. Ayo berfikir tentang pengujian ini didalam model regresi linear sederhana,

\[ Y_i = β_0 + β_1x_i + ϵ_i\].

Jika kita mengasumsikan hipotesis nolnya adalah benar, maka \(β_1 =0\) dan kita memiliki model,

\[Y_i = β_0 + ϵ_i\].

Dimodel ini, responnya tidak berdasarkan prediktornya. Jadi kita dapat berfikit bahws uji ini mengikuti jalan,

Dibawah \(H_0\) tidak terdapat hubungan linier yang signifikan diantara \(x\) dan \(y\).

Dibawah \(H_1\) terdapat hubungan linear yang signifikan diantara \(x\) dan \(y\).

Untuk contoh cars,

Dibawah \(H_0\) tidak terdapat hubungan linier yang signifikan diantara \(x\) dan \(y\).

Dibawah \(H_1\) terdapat hubungan linear yang signifikan diantara kecepatan dan jarak berhenti.

Lagi, uji itu dilihat pada keluaran dari summary(),

\[p-value = 1.4898365 × 10^{-12}\].

Dengan nilai p yang sangat rendah, kita menolak hipotesis nol tingkatan \(α\) dengan alasan apapun, katakanlah untuk contoh, \(α = 0.01\). Jadi kita mengatakan bahwa ada hubungan linear yang signifikan diantara kecepatan dan jarak berhenti. Perhatikan bahwa kita menekankan linearnya.

knitr::include_graphics("D:/Sem5/Computation/UTS/gambar1.jpg")

Di plot ini yang mensimulasikan data, kita dapat melihat secara jelas hubungan diantara \(x\) dan \(y\), lagipula ini bukanlah hubungan linear. Jika kita mencocokkan garis terhadap data ini, akan menjadi sangat datar. Hasil dari pengujian untuk \(H_0 : β_1 = 0\) melawan \(H_1 : β_1 ≠ 0\) memberikan nilai p yang besar, dikasus ini sebesae 0,7564548, jadi kita akan gagal untuk menolak dan mengatakan bahwa disana tidak ada hubungan linear yang signifikan diantara \(x\) dan \(y\). Nanti kita akan melihat bagaimana caranya untuk mencocokan kurva terhadap data ini menggunakan model “linear”, tapi untuk sekarang, sadari bahwa uji \(H_0 : β_1 = 0\) melawan \(H_1 : β_1 ≠ 0\) hanya dapat mendeteksi hubungan garis lurusnya saja.

1.6.3 Interval Keyakinan di R

Dengan menggunakan R, kita bisa mencari interval kepercayaan untuk \(β_0\) dan \(β_1\) dengan sangat mudah.

confint(stop_dist_model, level = 0.99)

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

Ini otomatis menghitung dengan interval kepercayaan sebesar 99% untuk kedua \(β_0\) dan \(β_1\), baris pertama untuk \(β_0\), baris kedua untuk \(β_1\).

Untuk contoh cars disaat menafsirkan interval ini, kita katakan, kita percaya bahwa 99% untuk kecepatan 1 mil per jam yang meningkat, peningkatan rata-rata jarak berhenti adalah diantara 2.8179187 dan 5.0468988 kaki, dimana merupakan interval untuk \(β_1\).

Perhatikan bahwa 99% interval kepercayaan ini tidak mengandung nilai hipotesis 0. Sejak itu tidak mengandung 0, ini setara dengan menolak uji \(H_0 : β_1 = 0\) melawan \(H_1 : β_1 ≠ 0\) di \(α = 0.01\), dimana kita sudah melihat itu sebelumnya.

Kamu harus sedikit curiga kepada interval kepercayaan untuk \(β_0\), sejak itu mencakup nilai negatif, yang sesuai dengan jarak berhenti negatif. Secara teknis, interpretasinya adalah kita percaya bahwa 99% jarak berhenti pada suatu mobil yang menempuh jarak 0 mil perjam berada diantara -35.7066103 dan 0.5484205 kaki, tetapi kita tidak sangat mempercayai itu, sejak kita sangat yakin bahwa itu akan menjadi non-negatif.

Catatan, kita dapat meng-ekstrak nilai spesifik dari keluaran ini dengan beberapa cara. Kode ini tidak berjalan, dan sebagai gantinya, kamu harus mengecek bagaimana itu bisa berhubungan dengan keluaran 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

Kita juga dapat memastikan bahwa perhitungan dari R telah menunjukkan interval untuk \(β_1\)

# 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 untuk Respon Rata-rata

Selain interval kepercayaan untuk \(β_0\) dan \(β_1\), terdapat dua estimasi interval umum lainnya yang digunakan dengan regresi. Yang pertama dipanggil confidence interval for the mean response. Seringkali, kita akan memperkirakan interval untuk rata-rata, \(E[Y | X = x]\) untuk nilai tertentu dari \(x\).

Di situasi ini kita menggunakan \(\hat y(x)\) sebagai estimasi kita dari \(E[Y | X = x]\). Kita memodifikasikan sedikit notasi kita untuk membuatnya jelas bahwa nilai prediksi adalah fungsi dari nilai \(x\).

\[\hat y(x) = \hat β_0 + \hat β_1 x\]

Ulangi lagi bahwa,

\[E[Y | X = x] = β_0 + β_1 x.\]

Jadi, \(\hat y(x)\) merupakan estimasi yang bagus sejak itu tidak bias:

\[E[\hat y(x)] = β_0 + β_1 x.\]

Kita bisa mendapatkan,

\[Var[\hat y(x)] = σ^2 (\frac {1}{n} = \frac {(x-x)^2}{S{xx}}).\]

Seperti estimasi lain yang telah kita lihat, \(\hat y(x)\) juga mengikuti distribusi normal. Sejak \(β_0\) dan \(β_1\) merupakan kombinasi linear dari variabel acak normal, begitupula dengan \(\hat y(x)\).

\[\hat y(x) ∼ N (β_0 + β_1x,σ^2(\frac{1}{n} + \frac {(x-x)^2}{S{xx}})).\]

Dan terakhir, sejak kita butuh untuk mengestimasi varians ini, kita sampai pada pada standar kesalahan dari estimasi kita,

\[SE[\hat y(x)] = s_e \sqrt {\frac {1}{n} + \frac {(x-\bar x)}{S{xx}}}\]

Lalu kita dapat menggunakan ini untuk mencari interval kepercayaan untuk respon rata-ratanya,

\[\hat y(x) ± t_{α/2,n-2} ⋅ s_e \sqrt {\frac {1}{n} + \frac {(x-\bar x)^2}{S{xx}} }\]

Untuk mencari tau interval kepercayaan untuk respon rata-rata menggunakan R, kita menggunakan fungsi predict(). Kita memberikan fungsi tersebut model kecocokkan kita sebagai data baru, menyimpannya sebagai bingkai data. (Ini penting, maka R diketahui sebagai nama dari variabel prediktor.) Disini, kita menemukan interva lkepercayaan untuk rata-rata jarak berhenti disaat sebuah mobil melaju 5 mil perjam dan disaat sebuah mobil melaju 21 mil perjam.

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 [ϵ]\] \[= σ^2 (\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}}) + σ^2 \] \[ = σ^2 (1 +\frac{1}{n}+\frac {(x-\bar x)^2}{S{xx}})\] \[\hat y(x) + ϵ ∼ N (β_0 + β_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_{α/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 seksi 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 = ∑_{i-1}^n (y_i - \bar y_i)^2 + ∑_{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("D:/Sem5/Computation/UTS/gambar2.jpg")

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 {∑_{i=1}^n ( y_i - \bar y)^2/1}{∑_{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 : β_1 = 0\] \[ Y_i = β_0 + ϵ_i\]

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

\[H_0 : β_1 ≠ 0\] \[ Y_i = β_1x_i + ϵ_i\]

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

\[F = \frac {∑_{i=1}^n (\hat y_i - \bar y)^2/1}{∑_{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 = β_0 + ϵ_i\) kepada data cars. Perhatikan bahwa \(\hat y = \bar y\) disaat \(Y_i = β_0 + ϵ_i\).

Pernyataan model lm(dist ~ speed, data = cars) menerapkan model \(Y_i = β_0 + β_1x_i + ϵ_i\).

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