Midterm Exam

1 Bab 8 Dugaan dalam Regresi Linier

Setelah membaca bab ini, diharapkan anda akan memiliki kemampuan untuk:

• Memahami distribusi perkiraan regresi.
• Membuat perkiraan interval untuk parameter regresi, respon rata-rata, dan prediksi.
• Uji signifikansi regresi.

Pada bagian akhir, kita akan mendefinisikan model regresi linier sederhana,

\[Y_i=\beta_0+\beta_1x_1+\epsilon_i\]

dimana \(\epsilon_i\backsim N(0,\sigma^2)\). Kemudian kita menggunakan pengamatan \((x_i,y_i)\), untuk \(i=1,2,...,n\), untuk mencari nilai \(\beta_0\) dan \(\beta_1\) yang diminimalkan

\[f(\beta_0,\beta_1)=\sum^n_{i=1}(y_i-(\beta_0+\beta_1x_i))^2\]

Kita memiliki nilai-nilai \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) yang telah kita temukan sebagai berikut:

\[\hat{\beta_1}=\frac{S_{xy}}{S_{xx}}=\frac{\sum^n_{i=1}(x_i-\bar{x})(y_i-\bar{y})}{\sum^n_{i=1}(x_i-\bar{x})^2}\]

\[\hat{\beta_0}=\bar{y}-\hat{\beta_1}\bar{x}\]

Kita juga memperkirakan \(\sigma^2\) menggunakan \(s^2_e\). Artinya, kita menemukan bahwa \(s_e\) ini adalah perkiraan dari \(\sigma\), dimana:

\[s_e=RSE=\sqrt{\frac{1}{n-2}\sum^n_{i=1}e^2_i}\]

yang mana RSE merupakan singkatan dari “Residual Standard Error” (Kesalahan Standar Sisa).

Ketika diterapkan pada data cars, kita memperoleh hasil sebagai 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

Pada bagian akhir, kita hanya membahas nilai-nilai Estimate, Residual Standard Error, dan Multiple R-squared. Pada bagian ini, kita akan membahas semua informasi mengenai Coefficients serta 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")

Untuk memulai, kita akan mencatat bahwa ada ekspresi lain yang setara untuk \(S_{xy}\) yang tidak kita lihat pada bagian akhir,

\[S_{xy}=\sum^n_{i=1}(x_i-\bar{x})(y_i-\bar{y})=\sum^n_{i=1}(x_i-\bar{x}y_i)\]

Ini mungkin kesetaraan yang mengejutkan (Cobalah untuk membuktikannya). Namun, hal tersebut akan sangat membantu untuk menggambarkan konsep pada bagian ini.

Perhatikan bahwa, \(\hat{\beta_1}\) adalah statistik sampel jika dihitung dengan data pengamatan yang ditulis di atas, seperti \(\hat{\beta_0}\).

Tetapi, dalam bagian ini akan lebih mudah untuk menggunakan \(\hat{\beta_1}\) dan \(\hat{\beta_0}\) sebagai variabel acak, akan tetapi, kita belum mengamati nilai untuk setiap \(Y_i\). Jika demikian, kita akan menggunakan notasi yang sedikit berbeda, menggantikan \(Y_i\) kapital dengan \(y_i\) huruf kecil.

\(\hat{\beta_1}=\frac{\sum^n_{i=1}(x_i-\bar{x})Y_i}{\sum^n_{i=1}(x_i-\bar{x})^2}\)

\(\hat{\beta_0}=\bar{Y}-\hat{\beta_1}\bar{x}\)

Pada bagian akhir kita berpendapat bahwa perkiraan parameter model yang tidak diketahui ini (\(\beta_0\) dan \(\beta_1\)) adalah baik karena kita memperolehnya dengan meminimalkan kesalahan. Sekarang kita akan membahas teorema Gauss-Markov yang membawa gagasan ini lebih lanjut, menunjukkan bahwa perkiraan ini sebanarnya adalah perkiraan yang “terbaik”, dari sudut pandang tertentu.

1.1 Teorema Gauss-Markov

Teorema Gauss-Markov memberitahu kita bahwa ketika memperkirakan parameter model regresi linier sederhana \(\beta_0\), dan \(\beta_1\), \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) yang kita dapatkan adalah perkiraan tidak bias linier terbaik. (Kondisi sebenarnya teorema Gauss-Markov lebih santai daripada model SLR)

Sekarang kita akan membahas linier, tidak bias, dan terbaik yang berhubungan dengan perkiraan ini.

1.1.1 Linier

Mengingat kembali, dalam aturan SLR, nilai-nilai \(x_i\) dianggap tetap dan jumlahnya diketahui. Kemudian, perkiraan linier adalah salah satu yang dapat ditulis sebagai kombinasi linier dari \(Y_i\). Dalam kasus \(\hat{\beta_1}\) kita dapatkan

\[\hat{\beta_1}=\frac{\sum^n_{i=1}(x_i-\bar{x})Y_i}{\sum^n_{i=1}(x_i-\bar{x})^2}=\sum^n_{i=1}k_iY_i=k_1Y_1+k_2Y_2+...+k_nY_n\]

dimana \(k_i=\frac{(x_i-\bar{x})}{\sum^n_{i=1}(x_i-\bar{x})^2}\).

Dengan cara yang sama, kita dapat menunjukkan bahwa \(\hat{\beta_0}\) dapat ditulis sebagai kombinasi linier dari \(Y_i\). Jadi, \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) keduanya adalah estimator linier.

1.1.2 Unbias

Sekarang kita tahu bahwa perkiraan kita adalah linier, seberapa baguskah perkiraan ini? Salah satu ukuran “goodness” (kebaikan) dari perkiraan ini adalah bias. Secara khusus, kita lebih suka perkiraan yang tidak bias, yang artinya nilai yang diharapkan adalah parameter yang sedang diperkirakan.

Pada kasus perkiraan regresi, kita memiliki

\[E[\hat{\beta_0}]=\beta_0\] \[E[\hat{\beta_1}]=\beta_1\]

Ini artinya, ketika menemukan kondisi model SLR, rata-rata perkiraan kita akan benar. Tetapi, seperti yang kita lihat pada bagian terakhir ketika menyimulasikan dari model SLR, itu bukan berarti bahwa setiap perkiraan individual akan benar. Hanya saja, jika kita mengulangi prosesnya sebanyak tak terhingga, rata-rata perkiraan akan benar.

1.1.3 The Best

Sekarang, jika kita membatasi diri terhadap perkiraan linier dan tidak bias, bagaimana cara kita menentukan perkiraan terbaik? Perkiraan dengan variansi minimum.

Pertama, perhatikan bahwa sangat mudah untuk membuat perkiraan untuk \(\beta_1\) yang memiliki variansi sangat rendah, tetapi tidak bias. Sebagai contoh, tentukan:

\[\hat{\theta}_{BAD}=5\].

Lalu, karena \(\hat{\theta}_{BAD}\) adalah nilai konstan, maka

\[Var[\hat{\theta}_{BAD}]=0\].

Namun karena,

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

kita mengatakan bahwa \(\hat{\theta}_{BAD}\) adalah estimator bias kecuali \(\beta_1=5\) yang tidak akan kita ketahui sebelumnya. Ini adalah perkiraan yang buruk (kecuai secara kebetulan \(\beta_1=5\)) meskipun memiliki kemungkinan variansi terkecil. Ini adalah bagian dari alasan kita membatasi diri pada perkiraan tidak bias. Apa gunanya sebuah perkiraan, jika memperkirakan kuantitas yang salah?

Sekarang, pertanyaan wajarnya adalah, apa variansi dari \(\hat{\beta_0}\) dan \(\hat{\beta_1}\)? Yaitu,

\[Var[\hat{\beta_0}]=\sigma^2(\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}})\]

\[Var[\hat{\beta_1}]=\frac{\sigma^2}{S_{xx}}\]

Ini mengukur variabilitas dari perkiraan dikarenakan peluang acak selama pengambilan sampel. Apakah ini “yang terbaik”? Apakah variansi sekecil ini yang kemungkinan kita bisa dapatkan? Anda hanya perlu mempercayai apa yang kami katakan karena menunjukkan bahwa ini benar berada di luar cakupan kursus ini.

1.2 Distribusi Sampling

Sekarang kita telah “mendefinisikan ulang” perkiraan untuk \(\hat{\beta_0}\) dan \(\hat{\beta_1}\) sebagai variabel acak, kita dapat mendiskusikan distribusi samplingnya, yang merupakan distribusi ketika statistik dianggap sebagai variabel acak.

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

kemudian, gabungkan semua hal di atas, kita dapatkan distribusi \(\hat{\beta_0}\) dan \(\hat{\beta_1}\).

Untuk \(\hat{\beta_1}\) kita tuliskan, \(\hat{\beta_1}=\frac{S_{xy}}{S_{xx}}=\frac{\sum^n_{i=1}(x_i-\bar{x})Y_i}{\sum^n_{i=1}(x_i-\bar{x})^2}\backsim N(\beta_1,\frac{\sigma^2}{\sum^n_{i=1}(x_i-\bar{x})^2})\)

atau lebih ringkas,

\[\hat{\beta_1}\backsim N(\beta_1,\frac{\sigma^2}{S_{xx}})\].

Dan untuk \(\hat{\beta_0}\),

\[\hat{\beta_0}\backsim N(\beta_0,\sigma^2(\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}))\]

Sejauh ini, kita gagal membuktikan sejumlah hasilnya. Kita tidak akan mengerjakan turunan yang membosankan dari distribusi sampel ini, tetapi kita akan membenarkan hasilnya menggunakan simulasi.

Catatan untuk pembaca saat ini: Turunan dan pembuktian ini mungkin bisa ditambahkan ke lampiran untuk kemudian hari. Anda juga dapat menemukan hasil ini di hampir semua buku regresi linier sederhana. Di UIUC, hasil ini kemungkinan besar akan disajikan di STAT 424 dan STAT425. Namun, karena Anda tidak diminta untuk melakukan penurunan jenis ini pada kursus kali ini, maka untuk saat ini hasil tersebut bisa dihilangkan.

1.2.1 Menyimulasikan Distribusi Sampling

Untuk memeriksa hasil di atas, kita akan menyimulasikan sampel dengan ukuran \(n=100\) dari model berikut

\[Y_i=\beta_0+\beta_1x_i+\epsilon_i\]

dimana \(\epsilon_i∼N(0,\sigma^2)\)

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

Kemudian, berdasarkan keterangan di atas, kita harus menemukan

\[\hat{\beta_1}\backsim N(\beta_1,\frac{\sigma^2}{S_{xx}})\]

\[\hat{\beta_0}\backsim N(\beta_0,\sigma^2(\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}))\].

Pertama, kita sebelumnya perlu memutuskan berapa nilai \(x\) kita untuk simulasi ini, karena nilai \(x\) dalam SLR kuantitasnya dianggap diketahui. Pilihan nilai \(x\) tidak tetap atau selalu berubah. Di sini kita juga menetapkan seed untuk pengacakan, dan menghitung \(S_{xx}\) yang akan kita perlukan untuk selanjutnya.

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

Kita juga akan memperbaiki nilai parameter kita.

beta_0 = 3
beta_1 = 6
sigma  = 2

Dengan informasi ini, kita mengetahui bahwa distribusi sampling seharusnya:

(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}\backsim N(6, 0.1176238)\)

dan

\(\hat{\beta_0}\backsim N(3, 0.04)\)

Yaitu,

\(E[\hat{\beta_1}]=6\) \(Var[\hat{\beta_1}]=0.1176238\)

dan

\(E[\hat{\beta_0}]=3\) \(Var[\hat{\beta_0}]=0.04\)

Sekarang kita menyimulasikan data dari model ini sebanyak 10.000 kali. Perhatikan bahwa cara ini bukan merupakan cara R untuk melakukan simulasi. Kita melakukan simulasi dengan cara ini untuk mendapatkan kejelasan. Sebagai contoh, kita bisa menggunakan fungsi sim_slr dari bagian sebelumnya. Kita juga menyimpan variabel di environment global, bukan membuat data frame untuk setiap kumpulan data simulasi 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 saat kita menyimulasikan data, kita memperoleh nilai perkiraan koefisien. Variabel beta_0_hats dan beta_1_hats sekarang menyimpan 10.000 nilai simulasi dari masing-masing \(\hat{\beta_0}\) dan \(\hat{\beta_1}\).

Kita verifikasi distribusi \(\hat{\beta_1}\) terlebih dahulu.

mean(beta_1_hats) # rata-rata empiris

## [1] 6.001998

beta_1            # rata-rata sebenarnya

## [1] 6

var(beta_1_hats)  # variansi empiris

## [1] 0.11899

var_beta_1_hat    # variansi sebenarnya

## [1] 0.1176238

Kita dapat lihat bahwa rata-rata maupun varianse empiris dan rata-rata maupun variansi sebenarnya sangatlah mirip. Kita juga akan memeriksa distribusi empiris bersifat normal. Untuk melakukannya, kita memplot beta_1_hats ke histogram, dan tambahkan kurva untuk distribusi sebenarnya dari \(\hat{\beta_1}\). Kita gunakan prob = TRUE untuk meletakkan histogram pada skala yang sama dengan kurva normal.

# catat bahwa kita perlu menggunakan 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)

Kemudian, kita ulangi prosesnya untuk \(\hat{\beta_0}\)

mean(beta_0_hats) # rata-rata empiris

## [1] 3.001147

beta_0            # rata-rata sebenarnya

## [1] 3

var(beta_0_hats)  # variansi empiris

## [1] 0.04017924

var_beta_0_hat    # variansi sebenarnya

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

Pada pelajaran simulasi ini, kita hanya bisa menyimulasikan sampel dengan jumlah yang terbatas. Untuk benar-benar memeriksa hasil distribusi, kita perlu mengamati sampel dengan jumlah yang tak terbatas. Namun, plot berikut ini menjelaskan bahwa jika kita melanjutkan proses simulasi ini, maka hasil empirisnya semakin mendekati apa yang kita harapkan.

par(mar = c(5, 5, 1, 1)) # sesuaikan margin plot, jika tidak "topi" tidak akan ditampilkan
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)) # sesuaikan margin plot, jika tidak "topi" tidak akan ditampilkan
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 Error

Jadi, sekarang kita yakin terhadap hasil distribusi berikut,

\[\hat{\beta_0}\backsim N(\beta_0,\sigma^2(\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}))\] \[\hat{\beta_1}\backsim N(\beta_1,\frac{\sigma^2}{S_{xx}})\]

Dengan menstandarkan hasil ini, kita temukan bahwa

\[\frac{\hat{\beta_0}-\beta_0}{SD[\hat{\beta_0}]}\backsim N(0,1)\]

dan

\[\frac{\hat{\beta_1}-\beta_1}{SD[\hat{\beta_1}]}\backsim N(0,1)\]

dimana

\[SD[\hat{\beta_0}]=\sigma\sqrt{\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}}\]

dan

\[SD[\hat{\beta_1}]=\frac{\sigma}{\sqrt{S_{xx}}}\].

Karena kita tidak mengetahui \(\sigma\) dalam praktiknya, kita harus memperkirakannya menggunakan \(s_e\), yang kita masukkan ke ekspresi yang ada untuk standar deviasi perkiraan kita.

Kedua ekspresi baru ini disebut dengan kesalahan standar yang merupakan perkiraan standar deviasi dari distribusi sampling.

\[SE[\hat{\beta_0}]=s_e\sqrt{\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}}\] \[SE[\hat{\beta_1}]=\frac{s_e}{\sqrt{S_{xx}}}\]

Sekarang, jika kita bagi dengan kesalahan standar (bukan standar deviasi), kita akan memperoleh hasil berikut yang memungkinkan untuk membuat interval kepercayaan dan melakukan uji hipotesis.

\[\frac{\hat{\beta_0}-\beta_0}{SE[\hat{\beta_0}]}\backsim t_{n-2}\] \[\frac{\hat{\beta_1}-\beta_1}{SE[\hat{\beta_1}]}\backsim t_{n-2}\]

Untuk melihatnya, pertama-tama perhatikan bahwa,

\[\frac{RSS}{\sigma^2}=\frac{(n-2)s^2_e}{\sigma^2}\backsim \chi^2_{n-2}\]

Dan juga ingat kembali variabel \(T\) didefinisikan sebagai,

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

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

Kita tulis,

\(T\backsim t_d\)

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

Kemudian, kita menggunakan trik klasik “dikalikan dengan 1” dan beberapa penyusunan ulang untuk sampai kepada,

\(\frac{\hat{\beta_1}-\beta_1}{SE[\hat{\beta_1}]}=\frac{\hat{\beta_1}-\beta_1}{\frac{s_e}{\sqrt{S_{xx}}}}\)

\(=\frac{\hat{\beta_1}-\beta_1}{\frac{s_e}{\sqrt{S_{xx}}}}.\frac{\frac{\sigma}{\sqrt{S_{xx}}}}{\frac{\sigma}{\sqrt{S_{xx}}}}\)

\(=\frac{\hat{\beta_1}-\beta_1}{\frac{\sigma}{\sqrt{S_{xx}}}}.\frac{\frac{\sigma}{\sqrt{S_{xx}}}}{\frac{s_e}{\sqrt{S_{xx}}}}\)

\(=\frac{\hat{\beta_1}-\beta_1}{\frac{\sigma}{\sqrt{S_{xx}}}}/\sqrt{\frac{s^2_e}{\sigma^2}}\)

\(=\frac{\hat{\beta_1}-\beta_1}{SD[\hat{\beta_1}]}/\sqrt{\frac{\frac{(n-2)s^2_e}{\sigma^2}}{n-2}}\backsim \frac{Z}{\sqrt{\frac{\chi^2_{n-2}}{n-2}}}\backsim t_{n-2}\)

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

Ingatlah bahwa distribusi \(t\) mirip dengan normal standar, tetapi dengan ekor (tails) yang lebih padat. Ketika derjat kebebasannya meningkat, distribusi \(t\) menjadi lebih seperti normal standar. Di bawah ini kita plot distribusi normal standar dengan dua contoh dari distribusi \(t\) dengan derajat kebebasan yang berbeda. Perhatikan bagaimana distribusi \(t\) dengan derajat kebebasan yang lebih besar akan lebih mirip dengan kurva normal standar.

# menentukan grid dari nilai x
x = seq(-4, 4, length = 100)

# memplot kurva untuk normal standar
plot(x, dnorm(x), type = "l", lty = 1, lwd = 2,
     xlab = "x", ylab = "Density", main = "Normal vs t Distributions")
# menambah kurva untuk distribusi t
lines(x, dt(x, df = 1), lty = 3, lwd = 2, col = "darkorange")
lines(x, dt(x, df = 10), lty = 2, lwd = 2, col = "dodgerblue")

# menambah 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 Kepercayaan untuk Kemiringan dan Intersep

Ingat bahwa interval kepercayaan untuk rata-rata sering kali berbentuk:

\(EST±CRIT.SE\)

atau

\(EST±MARGIN\)

dimana EST merupakan perkiraan untuk parameter yang diinginkan, SE adalah kesalahan standar dari perkiraan, dan MARGIN = CRIT.SE.

Kemudian, kita dapat membuat interval kepercayaan untuk \(\beta_0\) dan \(\beta_1\) menggunakan

\(\hat{\beta_0}±t_{a/2,n-2}.SE[\hat{\beta_0}]\) \(\hat{\beta_0}±t_{a/2,n-2}.s_e\sqrt{\frac{1}{n}+\frac{\bar{x}^2}{S_{xx}}}\)

dan

\(\hat{\beta_1}±t_{a/2,n-2}.SE[\hat{\beta_1}]\) \(\hat{\beta_1}±t_{a/2,n-2}.\frac{s_e}{\sqrt{S_{xx}}}\)

dimana \(t_{a/2, n-2}\) adalah nilai kritis seperti \(P(t_{n-2}>t_{a/2, n-2)}=\frac{a}{2}\).

1.5 Pengujian Hipotesis

Ingat bahwa statistik uji \(\text{TS}\) untuk sarana pengujian sering kali berbentuk:

\(\text{TS} = \frac{\text{EST} - \text {HYP}} {\text {SE}}\)

di mana \(\text {HYP}\) \(\text {EST}\) adalah perkiraan untuk parameter yang diminati, adalah nilai hipotesis dari parameter, dan \(\text{SE}\) adalah kesalahan standar perkiraan.

Jadi, untuk mengujinya

\(H_0: \beta_0 = \beta_ {00} \quad \text {vs} \quad H_1: \beta_0 \neq \beta_ {00}\)

kami menggunakan statistik uji

\(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, di bawah hipotesis nol, mengikuti a t distribusi dengan \(n - 2\) derajat kebebasan. Kita gunakan \(\beta_ {00}\) untuk menunjukkan nilai yang dihipotesiskan dari \(\beta_0.\)

Begitu pula untuk mengetes

\(H_0: \beta_1 = \beta_ {10} \quad \text {vs} \quad H_1: \beta_1 \neq \beta_ {10}\)

kami menggunakan statistik uji

\(t = \frac {\hat {\beta} _1- \beta_ {10}} {\text {SE} [\hat {\beta} _1]} = \frac {\hat {\beta} _1- \beta_ { 10}} {s_e / \sqrt {S_ {xx}}}\)

yang lagi-lagi, di bawah hipotesis nol, mengikuti distribusi t dengan \(n - 2\) derajat kebebasan. Sekarang kami menggunakan \(\beta_ {10}\) untuk menunjukkan nilai hipotesis dari \(\beta_1\) .

1.6 Contoh cars

Sekarang kita kembali ke cars contoh dari bagian terakhir untuk mengilustrasikan konsep ini. Kami pertama-tama menyesuaikan model menggunakan lm() kemudian menggunakan summary() untuk melihat hasil secara lebih rinci.

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 Pengujian dalam R

Sekarang kita akan membahas hasil yang ditampilkan yang disebut Coefficients. Pertama, ingatlah bahwa kita dapat mengekstrak 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

The names() Fungsi memberitahu kita apa informasi yang tersedia, dan kemudian kita menggunakan operator dan coefficient suntuk mengekstrak informasi yang kita tertarik. Dua nilai di sini harus segera akrab.

\[\hat {\beta} _0 = -17.5790949\]

dan

\[\hat {\beta} _1 = 3.9324088\]

yang merupakan perkiraan kami untuk parameter model \(\beta_0\) dan \(\beta_1\) .

Sekarang mari kita fokus pada keluaran baris kedua, yang relevan dengan \(\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

Sekali lagi, nilai pertama Estimate adalah

\[\hat {\beta} _1 = 3.9324088.\]

Nilai kedua, Std. Error, adalah kesalahan standar dari \(\hat {\beta} _1\),

\(\text {SE} [\hat {\beta} _1] = \frac {s_e} {\sqrt {S_ {xx}}} = 0,4155128.\)

Nilai ketiga, t-value, adalah nilai statistik uji untuk pengujian \(H_0: \beta_1 = 0 \quad {vs} \quad H_1: \beta_1 \neq 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|)\). beri kami nilai p dari pengujian itu.

\(\text {p-value} = 1,4898365 \times 10 ^ {- 12}\)

Perhatikan di sini, kami secara khusus menguji apakah \(\beta_1 = 0\) atau tidak .

Baris pertama keluaran melaporkan nilai yang 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

Singkatnya, kode berikut menyimpan informasi summary(stop_dist_model)coefficients dalam variabel baru (stop_dist_model_test_info), lalu mengekstrak setiap elemen menjadi variabel baru yang menjelaskan 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|)

Kemudian kita dapat memverifikasi beberapa ekspresi yang setara: statistik uji t untuk \(\hat {\beta} _1\) dan nilai p dua sisi yang terkait dengan statistik uji 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 Signifikansi Regresi, Uji-t

Kami berhenti sejenak untuk membahas pentingnya uji regresi . Pertama, perhatikan bahwa berdasarkan hasil distribusi di atas, kami dapat menguji \(\beta_0\) dan \(\beta_1\) terhadap nilai tertentu, dan melakukan pengujian satu dan dua sisi.

Namun, satu tes yang sangat spesifik,

\(H_0: \beta_1 = 0 \quad \text {vs} \quad H_1: \beta_1 \neq 0\)

paling sering digunakan. Mari kita pikirkan pengujian ini dalam kaitannya dengan model regresi linier sederhana,

\(Y_i = \beta_0 + \beta_1 x_i + \epsilon_i.\)

Jika kami mengasumsikan hipotesis nol benar, maka \(\beta_1 = 0\) dan kami memiliki model,

\(Y_i = \beta_0 + \epsilon_i.\)

Dalam model ini, respon tidak bergantung pada prediktor. Jadi kita bisa memikirkan tes ini dengan cara berikut,

• Di bawah \(H_0\) tidak ada hubungan linier yang signifikan antara \(x\) dan \(y\) .
• Di bawah \(H_1\) ada hubungan linier yang signifikan antara \(x\) dan \(y\) .

Untuk contoh cars,

• Di bawah \(H_0\) tidak ada hubungan linier yang signifikan antara kecepatan dan jarak berhenti.
• Di bawah \(H_1\) ada hubungan linier yang signifikan antara kecepatan dan jarak berhenti.

Sekali lagi, pengujian itu terlihat pada keluaran dari summary(),

\(\text {nilai-p} = 1,4898365 \times 10 ^ {- 12}.\)

Dengan nilai p yang sangat rendah ini, kami akan menolak hipotesis nol pada tingkat \(\alpha\) yang masuk akal , misalnya \(\alpha = 0,01\) . Jadi kami katakan ada hubungan linier yang signifikan antara kecepatan dan jarak berhenti. Perhatikan bahwa kami menekankan linier .

Dalam plot data simulasi ini, kita melihat hubungan yang jelas antara \(x\) dan \(y\) , namun ini bukan hubungan linier. Jika kita memasukkan garis ke data ini, itu sangat datar. Hasil pengujian untuk \(H_0: \beta_1 = 0 \quad \text{vs} \quad H_1: \beta_1 \neq 0\) memberikan nilai p yang besar, dalam hal ini \(0,7564548\) , jadi kami gagal menolak dan mengatakan bahwa tidak ada hubungan linier yang signifikan antara \(x\) dan \(y\) . Kita akan melihat nanti bagaimana menyesuaikan kurva ke data ini menggunakan model “linier”, tetapi untuk saat ini, sadari bahwa pengujian \(H_0: \beta_1 = 0\quad \text{vs} \quad H_1: \beta_1 \neq 0\) hanya dapat mendeteksi hubungan garis lurus.

1.6.3 Interval Keyakinan dalam R

Dengan menggunakan Rkita bisa dengan mudah mendapatkan interval kepercayaan untuk \(\beta_0\) dan \(\beta_1\) .

confint(stop_dist_model, level = 0.99)

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

Ini secara otomatis menghitung 99% interval keyakinan untuk \(\beta_0\) dan \(\beta_1\) , baris pertama untuk \(\beta_0\) , baris kedua untuk \(\beta_1\) .

Sebagai carscontoh saat menafsirkan interval ini, kami katakan, kami 99% yakin bahwa untuk peningkatan kecepatan 1 mil per jam, peningkatan rata-rata jarak berhenti adalah antara \(2.8179187\) dan \(5.0468988\) kaki, yang merupakan interval untuk \(\beta_1\) .

Perhatikan bahwa interval keyakinan 99% ini tidak mengandung nilai hipotesis 0. Karena tidak mengandung 0, ini setara dengan menolak pengujian \(H_0: \beta_1 = 0\) vs \(H_1: \beta_1 \neq 0\) di \(\alpha = 0.01\) , yang telah kami lihat sebelumnya.

Anda harus mencurigai interval keyakinan untuk \(\beta_0\) , karena mencakup nilai negatif, yang sesuai dengan jarak henti 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, kita dapat mengekstrak nilai tertentu dari keluaran ini dengan beberapa cara. Kode ini tidak dijalankan, dan sebagai gantinya, Anda harus memeriksa hubungannya dengan keluaran kode di atas.

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

Kami juga dapat memverifikasi bahwa kalkulasi yang bekerja Runtuk interval \(\beta_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 Keyakinan untuk Respon Rata-rata

Selain interval keyakinan untuk \(\beta_0\) dan \(\beta_1\) , ada dua perkiraan interval umum lainnya yang digunakan dengan regresi. Yang pertama disebut interval kepercayaan untuk respons rata-rata . Seringkali, kita menginginkan estimasi interval untuk mean, \(E [Y \mid X = x]\) untuk nilai \(x\) tertentu .

Dalam situasi ini kita menggunakan \(\hat {y} (x)\) sebagai perkiraan kita untuk \(E [Y \mid X = x]\) . Kami mengubah sedikit notasi kami untuk memperjelas bahwa nilai prediksi adalah fungsi dari nilai x .

\[\hat {y} (x) = \hat {\beta} _0 + \hat {\beta} _1 x\]

Ingat itu,

\[\text {E} [Y \mid X = x] = \beta_0 + \beta_1 x.\]

Jadi, \(\hat {y} (x)\) adalah perkiraan yang baik karena tidak bias:

\[\text {E} [\hat {y} (x)] = \beta_0 + \beta_1 x.\]

Kami kemudian bisa mendapatkan,

\[\text{Var}[\hat{y}(x)] = \sigma^2 \left(\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}\right).\]

Seperti perkiraan lain yang telah kita lihat, \(\hat {y} (x)\) juga mengikuti distribusi normal. Karena \(\hat {\beta} _0\) dan \(\hat {\beta} _1\) adalah kombinasi linier dari variabel acak normal, \(\hat {y} (x)\) juga.

\[\hat{y}(x) \sim N \left(\beta_0 + \beta_1 x, \sigma^2 \left(\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}\right) \right)\]

Dan terakhir, karena kita perlu memperkirakan varians ini, kita sampai pada kesalahan standar dari perkiraan kita,

\[\text {SE} [\hat {y} (x)] = s_e \sqrt {\frac {1} {n} + \frac {(x- \bar {x}) ^ 2} {S_ {xx}} }.\]

Kami kemudian dapat menggunakan ini untuk menemukan interval kepercayaan untuk respons rata-rata,

\[\hat {y} (x) \pm t _ {\alpha / 2, n - 2} \cdot s_e \sqrt {\frac {1} {n} + \frac {(x- \bar {x}) ^ 2 } {S_ {xx}}}\]

Untuk menemukan interval kepercayaan untuk respons rata-rata menggunakan R, kami menggunakan predict() fungsi. Kami memberikan fungsi model pas kami serta data baru, disimpan sebagai bingkai data. (Ini penting, agar R mengetahui nama variabel prediktor.) Di sini, kami menemukan interval kepercayaan untuk jarak berhenti rata-rata saat mobil melaju 5 mil per jam dan saat mobil melaju 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 Interval Prediksi untuk Pengamatan Baru

Terkadang kita menginginkan estimasi interval untuk pengamatan baru, \(Y\) , untuk nilai \(x\) tertentu . Ini sangat mirip dengan interval untuk respons rata-rata, \(\text {E} [Y \mid X = x]\) , tetapi berbeda dalam satu hal yang sangat penting.

Tebakan terbaik kami untuk observasi baru masih \(\hat {y} (x)\) . Perkiraan rata-rata masih merupakan prediksi terbaik yang bisa kita buat. Perbedaannya terletak pada jumlah variabilitas. Kita tahu bahwa pengamatan akan bervariasi tentang garis regresi yang sebenarnya menurut distribusi \(N (0, \sigma ^ 2)\) . Karena itu kami menambahkan faktor tambahan \(\sigma^2\) ke variabilitas perkiraan kami untuk menjelaskan variabilitas pengamatan tentang garis regresi.

\(\begin{aligned}\text{Var}[\hat{y}(x) + \epsilon] &= \text{Var}[\hat{y}(x)] + \text{Var}[\epsilon] \\[2ex] \sigma^2 \left(\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}\right) + \sigma^2 \\[2ex] \sigma^2 \left(1 + \frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}\right)\end{aligned}\)

\(\hat{y}(x) + \epsilon \sim N \left(\beta_0 + \beta_1 x, \ \sigma^2 \left(1 + \frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}\right) \right)\)

\(\text{SE}[\hat{y}(x) + \epsilon] = s_e \sqrt{1 + \frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}}\)

Kami kemudian dapat menemukan interval prediksi menggunakan,

\(\hat {y} (x) \pm t _ {\alpha / 2, n - 2} \cdot s_e \sqrt {1 + \frac {1} {n} + \frac {(x- \bar {x}) ^ 2} {S_ {xx}}}.\)

Untuk menghitung ini untuk satu set poin dalam Rpemberitahuan hanya ada sedikit perubahan dalam sintaks dari menemukan 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

Perhatikan juga bahwa kedua interval ini lebih lebar daripada interval kepercayaan yang sesuai untuk respons rata-rata.

1.9 Pita Keyakinan dan Prediksi

Seringkali kami ingin memplot kedua interval kepercayaan untuk respons rata-rata dan interval prediksi untuk semua kemungkinan nilai \(x\) . Kami menyebutnya band kepercayaan dan 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:

• Kami menggunakan ylim argumen untuk meregangkan sumbu - y dari plot, karena pita memperpanjang lebih jauh dari titik point.

• Kami menambahkan titik di titik \((\bar {x}, \bar {y})\) .

  • Ini adalah titik yang akan selalu dilewati garis regresi . (Pikirkan alasannya.)
  • Ini adalah titik di mana band keyakinan dan prediksi menjadi yang tersempit. Lihatlah kesalahan standar keduanya untuk memahami mengapa.

*Pita prediksi (titik-titik biru) kurang melengkung dibandingkan dengan pita kepercayaan (biru putus-putus). Ini adalah hasil dari faktor tambahan \(\sigma ^ 2\) yang ditambahkan ke varian pada nilai \(x\) apa pun .

1.10 Signifikansi Regresi, Uji-F

Dalam kasus regresi linier sederhana, uji \(t\) untuk signifikansi regresi sama dengan uji lain, uji \(F\) untuk signifikansi regresi. Kesetaraan ini hanya akan benar untuk regresi linier sederhana, dan di bagian selanjutnya kita hanya akan menggunakan uji \(F\) untuk signifikansi regresi.

Ingat dari bagian terakhir dekomposisi varians yang kita lihat sebelum menghitung \(R ^ 2\) ,

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

atau, singkatnya,

\[\text {SST} = \text {SSE} + \text {SSReg}.\]

Untuk mengembangkan uji \(F\), kita akan mengatur informasi ini ke dalam tabel ANOVA,

Sumber	Jumlah Kuadrat	Derajat Kebebasan	Rata-rata Kuadrat	F
Regresi	\(\sum^n_{i=1}(\hat{y}_i-y)^2\)	1	\(\frac{SSReg}{1}\)	\(\frac{MSReg}{MSE}\)
Eror	\(\sum^n_{i=1}(y_i-\hat{y}_i)^2\)	\(n-2\)	\(\frac{SSReg}{n-2}\)
Total	\(\sum^n_{i=1}(y_i-y)^2\)	\(n-1\)

ANOVA, atau Analisis Varians akan menjadi konsep yang sering kita bahas dalam kursus ini. Untuk saat ini, kita akan fokus pada hasil tabel, yaitu statistik \(F\) ,

\[F = \frac {\sum_ {i = 1} ^ {n} (\hat {y} _i - \bar {y}) ^ 2/1} {\sum_ {i = 1} ^ {n} (y_i - \hat {y} _i) ^ 2 / (n - 2)} \sim F_ {1, n - 2}\]

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

Ingat, dalam arti dari uji regresi, \(Y\) tidak tidak tergantung pada \(x\) di hipotesis nol.

\(H_0: \beta_1 = 0 \quad \quad Y_i = \beta_0 + \epsilon_i\)

Sedangkan pada hipotesis alternatif \(Y\) mungkin bergantung pada \(x\) .

\(H_1: \beta_1 \neq 0 \quad \quad Y_i = \beta_0 + \beta_1 x_i + \epsilon_i\)

Kita dapat menggunakan statistik F untuk melakukan tes ini.

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

Secara khusus, kami akan menolak nol ketika statistik \(F\) besar, yaitu, ketika ada kemungkinan rendah bahwa pengamatan bisa datang dari model nol secara kebetulan. Kami akan menghitung nilai p dengan menggunakan R

Untuk melakukan uji \(F\) di R Anda dapat melihat baris terakhir dari keluaran dari panggilan yang memberikan nilai statistik uji, derajat kebebasan yang relevan, serta nilai p dari pengujian. summary() yang disebut F-statistic

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 anova() fungsi tersebut untuk menampilkan informasi dalam 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 untuk pengujian. Anda harus memperhatikan bahwa nilai-p dari uji \(t\) adalah sama. Anda mungkin juga memperhatikan bahwa nilai statistik uji untuk uji \(t\) , \(9.46399\) , dapat dikuadratkan untuk mendapatkan nilai statistik \(F\) , \(89.5671065\) .

Perhatikan bahwa ada cara lain yang setara untuk melakukan ini R, di mana kita akan sering kembali 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\) ke data mobil. Perhatikan bahwa \(\hat {y} = \bar {y}\) ketika \(Y_i = \beta_0 + \epsilon_i\) .

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

Kami kemudian dapat menganggap penggunaan ini anova()sebagai perbandingan langsung kedua model. (Perhatikan bahwa kita mendapatkan nilai p yang sama lagi.)