# Chapter 8
Setelah membaca bagian ini, Anda akan dapat: - Memahami distribusi perkiraan regresi. - Membuat estimasi interval untuk parameter regresi, respons rata-rata, dan prediksi. - Uji signifikansi regresi. Bagian terakhir kami akan mendefinisikan model regresi linear sederhana, \[Yi=β0+β1xi+ϵi\] Dimana, \(ϵi∼N(0,σ2)\). Kemudian menggunakan pengamatan (xi,yi), for i=1,2,…n, untuk menemukan nilai β0 dan β1 which minimized \[f(β0,β1)=∑^n_{i=1}(yi−(β0+β1xi))^2\]. Kami sebut nilai ini \(\hat{\beta_0}\) and \(\hat{\beta_1}\) , dimana dapat ditemukan dengan. \[\hat{\beta_1}={S_{xy} \over S_{xx}} = {∑^n_{i=1}(x_{i}-\bar{x})(y_{i}-\bar{y}) \over ∑^n_{i=1}(x_{i}-\bar{x})^2}\] \[\hat{\beta_0}=\bar{y}-\hat{\beta_1}\bar{x}\].
Kami juga memperkirakan \(\sigma^2\) menggunakan \(s^2_e\). Dengan kata lain, kita menemukan \(s_e\) dalam memperkirakan \(\sigma\), dimana ; \[s_e = \text{RSE} = \sqrt{\frac{1}{n - 2}\sum_{i = 1}^n e_i^2}\] Dimana kita juga menyebut RSE sebagai “kesalahan standar residual.” Ketika diterapkan pada data \(cars\), kami memperoleh hasil berikut:
##
## 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 terakhir, kami hanya membahas estimasi, kesalahan standar residual, dan beberapa nilai R-kuadrat. Di bagian ini, kita akan membahas semua informasi Koefisien serta 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") Untuk memulai, kami akan mencatat bahwa ada ekspresi setara lain untuk \(S_{XY}\) yang kami tidak melihat di bagian terakhir, \[S_{xy}= \sum_{i = 1}^{n}(x_i - \bar{x})(y_i - \bar{y}) = \sum_{i = 1}^{n}(x_i - \bar{x}) y_i.\] Ini mungkin kesetaraan yang mengejutkan. (Mungkin mencoba untuk membuktikannya.) Namun, ini akan berguna untuk mengilustrasikan konsep di bagian ini. Perhatikan bahwa \(\hat{\beta_1}\) adalah sampel statistik jika dihitung dengan data yang diamati seperti yang tertulis diatas , sebagaiamana \(\hat{\beta_0}\). Namun dibagian ini akan sering menggunakan keduanya \(\hat{\beta_1}\) dan \(\hat{\beta_0}\) sebagai variabel random, artinya, kita belum mengamati nilai-nilai untuk masing-masing \(Y_i\). Ketika ini terjadi, kita akan menggunakan notasi yang sedikit berbeda, penggantian model \(Y_i\) untuk huruf kecil \(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}\] Bagian terakhir kami berpendapat bahwa perkiraan parameter model yang tidak diketahui ini \(\beta_0\) dan \(\beta_1\) karena kami mendapatkannya dengan meminimalkan kesalahan. Kami sekarang akan membahas THeorema Gauss-Markov yang mengambil ide ini lebih lanjut, yang menunjukkan bahwa perkiraan ini sebenarnya adalah perkiraan “terbaik”, dari sudut pandang tertentu.
Theorema Gauss-Markov memberi tahu kita bahwa ketika memperkirakan parameter model regresi linier sederhana \(\beta_{0}\) dan \(\beta_{1}\), \(\beta_{0}\) dan \(\beta_{1}\) yang kami turunkan adalah perkiraan linear tidak bias terbaik, atau BLUE singkatnya. (Kondisi aktual Gauss-Markov lebih sederhana daripada model SLR.) Kami sekarang akan membahas linear, tidak bias, dan terbaik seperti yang berkaitan dengan perkiraan ini.
Recall, dalam pengaturan SLR bahwa nilai \(X_i\) dianggap sebagai jumlah tetap dan diketahui. Kemudian perkiraan linear adalah salah satu yang dapat ditulis sebagai kombinasi linear dari \(Y_i\), dalam kasus \(\hat{\beta}_1\) kita lihat \[\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={(x_{i}- \bar{x}) \over \sum_{i = 1}^n (x_{i}- \bar{x})^2}\] Dengan cara yang sama, kita bisa menunjukkan bahwa \(\hat\beta_0\) dapat ditulis sebagai kombinasi linear \(Y_i\) dengan demikian \(\hat\beta_0\) dan \(\hat\beta_1\) adalah pengukur linear.
Sekarang kita tahu bahwa perkiraan kita linier, seberapa baik perkiraan ini? Salah satu ukuran “kebaikan” dari perkiraan adalah biasnya. Secara khusus, kita lebih suka perkiraan yang tidak bias, yang berarti nilai yang diharapkan adalah parameter yang diperkirakan.
Dalam kasus perkiraan regresi, kita memiliki,\[\begin{aligned} \text{E}[\hat{\beta}_0] &= \beta_0 \\ \text{E}[\hat{\beta}_1] &= \beta_1. \end{aligned}\] Ini memberi tahu kita bahwa, ketika kondisi model SLR terpenuhi, rata-rata perkiraan kita akan benar. Namun, seperti yang kita lihat bagian terakhir ketika mensimulasikan dari model SLR, itu tidak berarti bahwa setiap perkiraan individu akan benar. Hanya saja, jika kita mengulangi prosesnya dalam jumlah yang tak terbatas, rata-rata perkiraannya akan benar.
Sekarang, jika kita membatasi untuk kedua linear dan perkiraan tidak bias, bagaimana kita mendefinisikan perkiraan terbaik? Perkiraan dengan varians minimum. Pertama mencatat bahwa sangat mudah untuk membuat perkiraan untuk \(\beta_1\) yang memiliki varians sangat rendah, tetapi tidak bias. Misalnya, tentukan: \[\hatθ_{BAD}=5\] Lalu, sejak \(\hatθ_{BAD}\) nilai konstant, \[Var[\hatθ_{BAD}]=0\] Namun, \[E[\hatθ_{BAD}]=0\] kita katakan bahwa \(\hatθ_{BAD}\) adalah penaksir bias kecuali \(\beta_{1}=5\) yang tidak diketahui sebelumnya. Untuk alasan ini , ini merupakan persamaan yang buruk (kecuali secara kebetulan \(\beta_{1}=5\)) meskipun variansinya sekecil mungkin. Ini adalah bagian dari alasan membatasi pada perkiraan yang tidak bias. Apa yang baik adalah perkiraan, jika memperkirakan jumlah yang salah?
Jadi sekarang, pertanyaan nya adalah, apa varians dari \(\hat\beta_0\) dan \(\hat\beta_0\)? Mereka adalah, \[{Var[\hat\beta_0]= σ^2({1\over n}+{\bar{x}^2 \over S_{xx}})}\] \[Var[\hat\beta_{1}]={{σ^2}\over{S_xx}}\] Ini mengukur variabilitas perkiraan karena peluang acak selama pengambilan sampel. Apakah ini “yang terbaik”? Apakah varians ini sekecil yang kita bisa dapatkan? Anda hanya perlu mengambil kata-kata kami untuk itu bahwa mereka karena menunjukkan bahwa ini benar adalah di luar lingkup kursus ini.
Sekarang kita telah “mendefinisikan ulang” perkiraan untuk \(\hat\beta_0\) dan \(\hat\beta_1\) sebagai variabel acak, kita dapat mendiskusikan distribusi pengambilan sampel mereka, yang merupakan distribusi ketika statistik dianggap sebagai variabel acak.
Karena keduanya \(\hat\beta_0\) dan \(\hat\beta_1\) adalah kombinasi linear dari \(Y_i\) dan masing-masing \(Y_i\) adalah distribusi normal, maka keduanya \(\hat\beta_0\) dan \(\hat\beta_1\) juga mengikuti distribusi normal. Kemudian, menempatkan semua hal di atas bersama-sama, sehingga mendapatkan distribusi \(\hat\beta_0\) dan \(\hat\beta_1\). Untuk \(\hat\beta_1\) dikatakan. \[\hat{\beta_1}={S_{xy} \over S_{xx}} = {∑^n_{i=1}(x_{i}-\bar{x})Y_i \over ∑^n_{i=1}(x_{i}-\bar{x})^2}∼N ({\beta_1}, {{\sigma^2} \over {\sum_{i = 1}^{n}(x_{i}-\bar{x})^2 }})\] atau lebih ringkas, \[{\hat\beta_1} ∼N ({\hat\beta_1}, {{\sigma^2}\over {S_{xx}}})\] Dan untuk \(\hat\beta_0\) \[\hat\beta_0=\bar(Y)-\hat\beta_{1}\bar(x)∼N({\beta_0}, {{\sigma^2}\sum_{i = 1}^{n}x^2_1 \over {\sum_{i = 1}^{n}(x_{i}-\bar{x})^2 }})\] atau lebih ringkas, \[\hat\beta_0∼N(\beta_0, \sigma^2({1\over n}+{\bar{x}^2\over {S_{xx}}}))\] Pada titik ini kami telah mengabaikan untuk membuktikan sejumlah hasil ini. Alih-alih bekerja melalui turunan dari distribusi pengambilan sampel ini, kita justru akan membenarkan hasil ini menggunakan simulasi. Catatan untuk pembaca saat ini: Derivasi dan bukti ini dapat ditambahkan ke lampiran di lain waktu. Anda juga dapat menemukan hasil ini di hampir semua buku teks regresi linear standar. Di UIUC, hasil ini kemungkinan akan disajikan di STAT 424 dan STAT 425. Namun, karena Anda tidak akan diminta untuk melakukan derivasi jenis ini, bagian ini dihilangkan.
Verifikasi hasil diatas, kita akan melakukan simulasi sampel dengan ukuran \(n=100\) dari model \[Y_{1}={\beta_0}+{\beta_{1}x_{1}+ϵ_i}\] Dimana \(ϵ_i∼N(0,\sigma^2)\). Dalam kasus ini, parameter yang diketahui antara lain : - \(\beta_0=3\) - \(\beta_1=6\) - \(\sigma^2=4\) Lalu, berdasarkan di atas, kita dapat menemukan \[\hat\beta_1 ∼
N(\beta_1, {{\sigma^2}\over {S_{xx}}})\] \[{\hat\beta_0} ∼N (\beta_0,\sigma^2 {1\over {n}}+{\bar x^2\over S_{xx}})\] Pertama kita perlu memutuskan sebelumnya apa nilai \(x\) kita
yang akan disimulasikan, nilai dalam SLR juga dianggap sebagai jumlah yang diketahui. Pilihan nilai \(X\) bebas Di sini kita juga menetapkan benih untuk pengacakan, dan menghitung \(S_{XX}\) yang akan kita butuhkan ke depan.
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 menentukan nilai Parameter
Dengan informasi ini, kita tahu distribusi pengambilan sampel harus:
## [1] 0.1176238
## [1] 0.04
\[\hat\beta_1∼N(6,0.1176238)\] dan \[\hat\beta_0∼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.4\] sekarang kita akan mensimulasikan data 10.000 kali dari model ini . Perhatikan ini mungkin bukan cara R untuk melakukan simulasi. Kita melakukan simulasi dengan cara ini dalam upaya kejelasan. Misalnya, kita bisa menggunakan fungsi sim_slr() dari bagian sebelumnya. Kita juga hanya menyimpan variabel di lingkungan global untuk membuat bingkai data untuk setiap set data baru yang disimulasikan.
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 memperoleh nilai koefisien variabel yang diperkirakan. Variabel beta_0_hats dan beta_1_hats sekarang simpan 10.000 simulasi nilai masing-masing \(\hat\beta_0\) dan \(\hat\beta_1\) Pertama memverifikasi distribusi \(\hat\beta_1\)
## [1] 6.001998
## [1] 6
## [1] 0.11899
## [1] 0.1176238
Kita melihat bahwa cara dan varians empiris sangat mirip. Kita juga memverifikasi bahwa distribusi empiris adalah normal. Untuk melakukannya, kami merencanakan histogram beta_1_hats, dan menambahkan kurva untuk distribusi \(\hat\beta_1\). Menggunakan prob = TRUE untuk menempatkan histogram pada 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 mengulang proses untuk \(\hat\beta_0\)
## [1] 3.001147
## [1] 3
## [1] 0.04017924
## [1] 0.04
## [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) Dalam studi simulasi ini, hanya mensimulasikan sejumlah sampel terbatas. Untuk benar-benar memverifikasi hasil distribusi, kita perlu mengamati sejumlah sampel yang tak terbatas. Namun, plot berikut harus memperjelas bahwa jika terus mensimulasikan, hasil empiris akan semakin dekat dan lebih dekat dengan apa yang diharapkan.
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) ## 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"))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\].
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\).
carsSekarang menggunakan cars dari bagian sebelumnya untuk mengilustrasi konsepnya. Awalnya gunakan model lm() lalu gunakan summary() untuk melihat hasilnya dengan lebih detail.
##
## 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
Kita sekarang akan mendiskusikan tentang hasil yang ditampilkan yang dipanggil Coefficients. Pertama, mengulang lagi bahwa kita dapat meng-ekstrak informasi ini secara langsung.
## [1] "call" "terms" "residuals" "coefficients"
## [5] "aliased" "sigma" "df" "r.squared"
## [9] "adj.r.squared" "fstatistic" "cov.unscaled"
## 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\).
## 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\).
## 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.
## [1] 9.46399
## [1] 9.46399
## [1] 1.489836e-12
## [1] 1.489836e-12
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,
Untuk contoh cars,
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.
## Warning: package 'OpenImageR' was built under R version 4.0.2
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.
RDengan menggunakan R, selang kepercayaan \(\beta_0\) dan \(\beta_1\) sangatlah mudah untuk didapatkan.
## 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.
## 0.5 % 99.5 %
## -35.7066103 0.5484205
## [1] -35.70661
## [1] 0.5484205
## 0.5 % 99.5 %
## (Intercept) -35.70661 0.5484205
## 0.5 % 99.5 %
## 2.817919 5.046899
## [1] 2.817919
## [1] 5.046899
## 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
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
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.
## 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.
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.
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,
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.
##
## 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.
## 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.
## 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.)