Email           :
Instagram   : https://www.instagram.com/marvis.zerex/
RPubs          : https://rpubs.com/invokerarts/
Linkedin     : https://www.linkedin.com/in/jeffry-wijaya-087a191b5/
Majors         : Business Statistics
Address      : ARA Center, Matana University Tower
                         Jl. CBD Barat Kav, RT.1, Curug Sangereng, Kelapa Dua, Tangerang, Banten 15810.



1 Pengenalan

Seperti yang Saya sudah katakan sebelumnya, Model Regresi digunakan untuk menunjukan hubungan antar variabel dengan garis yang memnuhi data observasi. Model regresi linear menggunakan garis lurus, sedangkan model regresi logistik dan linear menggunakan garis melengkung. Regresi memungkinkan Anda untuk memperkirakan variabel tetap berubah seperti variabel bebasnya berubah.

Dalam bagian ini, kita akan fokus dalam bagaimana menggunakan Regresi Linear Sederhana untuk memperkirakan hubungan antara 2 variabel-variabel kuantitatif. Kita juga dapat mengugnakan regresi linear sederhana untuk mengetahui:

  • Sebetapa kuat hubungan antara kedua variabel.
  • Nilai dari variabel tetap pada suatu nilai variabel bebas.



Contoh:

Asumsikan Anda adalah peneliti sosial yang tertarik dengan hubungan antara pendapatan dan kesenangan pelanggan. Anda melakukan survei kepada 500 orang yang berpenghasilan dengan interval dari $15k hingga $75k dan bertanya kepada mereka tingkat kebahagiaan dari skala 1 sampai 10. Kedua variabel Anda baik variabel tetap (pendapatan) dan variabel bebas (kesenangan) bersifat kuantitatif, jadi anda dapat melakukan analisis regresi untuk melihat apakah ada hubungan linear diantara kedua variabel tersebut.


2 Model Umum

Sebuah model regresi linear menasumsikan bahwa ada hubungan linear diantara ekspetasi kondisi dari kedua variabel kuantitatif:

where,

  • \(Y\) biasa disebut juga prediktor, penjelas atau variabel tetap.
  • \(X\) biasa disebut juga respon, hasil atau variabel bebas.
  • \(\beta_1\) adalah tingkat intercept atau koefisien.
  • \(\beta_2\) adalah tingkat slope atau koefisien.
  • \(\hat{Y}\) adalah nilai estimasi dari \(Y\) yang dihasilkan oleh \(X\).
  • \(b_1\) dan \(b_2\) nilai estimasi koefisien.

3 Regresi Linear Sederhana dalam R

Data untuk contoh ini terdapat pada package R yaitu PoEdata. Pertama-tama, Anda harus menginstall package “PoEdata”, ketik perintah berikut pada bagian Console Rstudio Anda:

install.packages("devtools")                         # Jika belum di install
library(devtools)                                    # Meng-load packagenya
install_git("https://github.com/ccolonescu/PoEdata") # Mengintall `PoEdata`

Kemudian, pertimbangkan data set:

library(PoEdata)                                     # Meng-load package `PoEdata`
library(DT)                                          # Menggunakan fungsi datatable 
data("food", package="PoEdata")
datatable(food)

Akan selalu menjadi ide yang baik untuk memeriksa data secara visual menggunakan diagram scatter, dimana dapat dibuat menggunakan fungsi plot(). Gambar dibawah adalah diagram scatter dari pengeluaran makanan berdasarkan penghasilan, dapat dilihat bahwa ada hubungan positif antara pendapatan dan pengeluaran makanan.

plot(food$income, food$food_exp, 
     ylim=c(0, max(food$food_exp)),
     xlim=c(0, max(food$income)),
     xlab="weekly income in $100", 
     ylab="weekly food expenditure in $", 
     type = "p")

Tentang \(R^2\)


3.1 Memperkirakan Regresi Linear

Untuk data pengeluaran makanan, model regresinya akan menjadi

\[ \begin{align} \tag{1} Y&=\beta_1+\beta_2 *X+e \\ food_{exp}&=\beta_1+\beta_2*income+e \end{align} \]

Fungsi R untuk memperkirakan model regresi linear adalah lm(y~x, data) dimana, jika digunakan dengan sendirinya tidak menunjukan hasil yang jelas; Sangat berguna untuk memberikan nama tersebut model sepertimod1, kemudian menunjukan hasilnya menggunakan summary(mod1).

mod1 <- lm(food_exp ~ income, data = food)           # Model Linear
smod1<- summary(mod1)                                # Mendapatkan Summary 
smod1
## 
## Call:
## lm(formula = food_exp ~ income, data = food)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -223.025  -50.816   -6.324   67.879  212.044 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   83.416     43.410   1.922   0.0622 .  
## income        10.210      2.093   4.877 1.95e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 89.52 on 38 degrees of freedom
## Multiple R-squared:  0.385,  Adjusted R-squared:  0.3688 
## F-statistic: 23.79 on 1 and 38 DF,  p-value: 1.946e-05

Jika anda tertarik hanya dengan beberapa hasil dari regresi, seperti estimasi koefisien, Anda dapat mendapatkannya dengan menggunakan perintah khusus. seperti coef().

b1 <- coef(mod1)[[1]]                                # nilai estimasi dari  β1
b2 <- coef(mod1)[[2]]                                # nilai estimasi dari  β2
b1
## [1] 83.416
b2
## [1] 10.20964

Daru hasil tersebut, kita dapat menyatakan persamaan estimasi regresi linear kita menjadi:

\[ \begin{align} \tag{2} \hat{Y}&=b_1+b_2*X \\ \hat{food_{exp}}&=83.416+10.20964*income \end{align} \]

Tingkat intercept yaitu, \(\beta_1\) , biasanya hanya sedikit penting dalam model ekonometrik; kita akan sering lebih tertarik pada tingkat slope yaitu, \(\beta_2\). Nilai perkitaan dari \(\beta_2\) menunjukan bahwa pengeluaran makanan untuk kebanyakan keluarga meningkan sebesar 10.209643 ketika penghasilan keluarga naik sebesar 1, dimana dalam kasus ini adalah $100. Fungis R abline() menambahkan garis regresi pada diagram scatter yang sudah dibuat sebelumnya menjadi seperti dibawah ini:

plot(food$income, food$food_exp, 
     xlab="weekly income in $100", 
     ylab="weekly food expenditure in $", 
     type = "p")
abline(b1,b2)

Bagaimana cara mengambil berbagai kesimpulan regresi? Untuk mendapatkan sebuah kesimpukan, Anda cukup mengetik nama model Anda, diikuti dengan tanda $ dan nama dari hasil yang ingin anda cari. Sebagai contoh, juga kita ingin mengetahui vektor koefisien dari mod1, kita menulisnya sebagaimod1\$coefficients dan smod1\$coefficients:

mod1$coefficients
## (Intercept)      income 
##    83.41600    10.20964
smod1$coefficients
##             Estimate Std. Error  t value     Pr(>|t|)
## (Intercept) 83.41600  43.410163 1.921578 6.218242e-02
## income      10.20964   2.093264 4.877381 1.945862e-05
names(mod1)                                          
names(smod1)                                         

Seperti yang kita lihat sebelumnya, beberapa hasil dapat dicari menggunakan fungsi khusus seperti coef(mod1), resid(mod1), fitted(mod1), and vcov(mod1).


3.2 Perkiraan dengan Model Regresi Linear Sederhana

Perkiraan tingkat regresi, \(b_1\) dan \(b_2\) memungkinkan kita untuk memprediksi ekspetasi pengeluaran makanan untuk setiap penghasilan yang diketahui. Semua yang kita butuhkan hanya cukup memasukan estimasi nilai parameter dan penghasilan yang ada pada persamaan seperti persamaan ke-2. Sebagai contoh, nilai ekspetasi dari food_exp Untuk penghasilan $2000 dihitunga dalam persamaan ke-3 . (Jangan lupa untuk membagi nilai penghasilan dengan 100 karena data variabel yang kita gunakan dihitung dalam ratusan dolar)

\[ \begin{align} \tag{3} \hat{food_{exp}}&=83.416+10.20964*20\\ & = $ 287.608861 \end{align} \]

Namun, R akan melakukan kalkulasi dengan fungsi predict(). Mari kita coba memperluas contoh yang kita punya dengan menggunakan lebih dari satu penghasilan untuk memrediksi pengeluaran makanan, anggap penghasilan = $2000, $2500, dan $2700. Fungsi predict() dalam R membutuhkan nilai baru dari variabel bebas dibentuk dalam suatu data frame. Bahkan ketika kita hanya ingin memprediksi 1 penghasilan, kita membutuhkan struktur data framte yang sama. Dalam R, sebuah kelompok angka dinyanakan dalam struktur c(). Perintah berikut dapat mempermudah anda mengerti.

mod1 <- lm(food_exp~income, data=food)
newx <- data.frame(income = c(20, 25, 27))
yhat <- predict(mod1, newx)
yhat                                                 # menunjukan hasil
##        1        2        3 
## 287.6089 338.6571 359.0764

3.3 Nilai Koefisien Regresi

Koefisien regresi \(b_1\) dan \(b_2\) adalah variabel random, karena mereka bergantung kepada sampel. Dalam kasus ini, kita mengulang sampel pada nilai koefisien regresi. Kita akan membuat angka dari subsample random dari data food dan menghitung ulang \(b_1\) and \(b_2\). Subsample random dapat dibuat menggunakan fungsi sample(), contoh begitu mengilustrasikan hanya untuk mencari \(b_2\).

N <- nrow(food)   # mengetahui jumlah data observasi dalam dataset
C <- 50           # jumlah subsample yang diharapkan
S <- 38           # jumlah sample yang diharapkan

sumb2 <- 0
for (i in 1:C){   # pengulangan jumlah dari subsample
  set.seed(3*i)   # kode yang berbeda untuk setiap subsample
  subsample <- food[sample(1:N, size=S, replace=TRUE), ]
  mod2 <- lm(food_exp~income, data=subsample)
  #jumlah b2 untuk semua subsample:
  sumb2 <- sumb2 + coef(mod2)[[2]]
}
print(sumb2/C, digits = 3)
## [1] 9.89

Hasilnya, \(b_2=9.89\), adalah rata-rata dari estimasi 50 \(b_2\).


3.4 Estimasi Varians and Kovarians

Banyak proses yang diperlukan untuk mengestimasi varians dan kovarians dari sebuah koefisien regresi. R menyatukannya dengan fungsi vcov():

(varb1 <- vcov(mod1)[1, 1])
## [1] 1884.442
(varb2 <- vcov(mod1)[2, 2])
## [1] 4.381752
(covb1b2 <- vcov(mod1)[1,2])
## [1] -85.90316

3.5 Perbedaan Varians dan Kovarians Tinggi Maupun Rendah

  • Varians
    • Tinggi : Semakin tinggi nilai varians berarti semakin dekat data-data dengan rata-ratanya
    • Rendah : Semakin rendah nilai varians berarti semakin jauh data-data dengan rata-ratanya
  • Kovarians
    • Tinggi : Semakin tinggi nilai kovarians berarti menunjukan bahwa ada hubungan yang kuat antar variabel
    • Rendah : Semakin rendah nilai kovarians berarti menunjukan bahwa ada hubungan yang lemah antar variabel

4 Hubungan Non-Linear

Terkadang, diagram scatter atau beberapa perkiraan secara teori menunjukan hubungan non-linear. Hubungan non-linear yang paling sering ditemukan ada 2, yaitu yang berhubungan dengan logaritma dari variabel tetap atau bebas dan persamaan polinomial. Model kuadrat membutuhkan kuadrat dari variabel bebas.

\[ \begin{align} \tag{4} Y=\beta_1+\beta_2*X^2+e\\ \end{align} \] Dalam R, variabel bebas melibatkan operasi matematika yang meliputi fungsi regresi dengan fungsi I(). Contoh berukut menggunakan dataset br dari packkage PoEdata, yang mana mencakup harga penjualan dan atribut dari 1080 rumah di Baton Rouge, LA. Harga dalam penjualan dihitung dalam dolar, dan dqft merupakan luas tanah dalam satuan kaki kuadrat.

data(br)
mod2 <- lm(price~I(sqft^2), data=br)
b1 <- coef(mod2)[[1]]
b2 <- coef(mod2)[[2]]
sqftx=c(2000, 4000, 6000)            # nilai yang diberikan untuk luas tanah
pricex=b1+b2*sqftx^2                 # harga yang sesuai dengan luas tanah yang diberikan 
DpriceDsqft <- 2*b2*sqftx            # efek marginal dalam luas tanah pada harga
elasticity=DpriceDsqft*sqftx/pricex 
b1; b2; DpriceDsqft; elasticity      # menunjukan hasil
## [1] 55776.57
## [1] 0.0154213
## [1]  61.68521 123.37041 185.05562
## [1] 1.050303 1.631251 1.817408

Sekarang kita ingin menunjukan diagram scatter dan melihat bagaimana persamaan kuadrat mencocokan data. R Chuck selanjutnya menunjukan dua cara untuk membuat grafik tersebut. Pertama-tama, cukup menunjukan persamaan kuadrat pada diagram scatter, menggunakan fungsi curve(); kedua, gunakan garis persamaan, yang mana dibutuhkan pengurtan data secata menaik berdasarkan sqft sebelum model regresi dievaluasi

plot(br$sqft, br$price, xlab="Total square feet", 
     ylab="Sale price, $", col="grey")
#menambahkan persamaan kuadrat pada scatterplot:
curve(b1+b2*x^2, col="red", add=TRUE) 

Cara alternatif untuk menggambar kurva yang cocok adalah:

ordat <- br[order(br$sqft), ] #mengurutkan data berdasarkan `sqft`
mod2 <- lm(price~I(sqft^2), data=ordat)
plot(br$sqft, br$price, 
     main="Dataset ordered after 'sqft' ", 
     xlab="Total square feet", 
     ylab="Sale price, $", col="grey")
lines(fitted(mod31)~ordat$sqft, col="red")

Model log-linear meregresi log dari variabel tetap pada hasil linear dari variabel bebas. (kecuali sudah ditentukan, notasi log menunjuk kepada logaritma natural, mengikuti konvensi umum pada bidang ekonomi):

\[ \begin{align} \tag{5} Log(Y)=\beta_1+\beta_2*X+e\\ \end{align} \]

Salah satu alasan menggunakan logaritma pada variabel bebas adalah membuat distribusinya mendekat kepada distribusi normal. Mari kita tunjukan histogram dari harga dan log(harga) untuk membandingkannya

hist(br$price, col='grey')

hist(log(br$price), col='grey')

Perbedaan setelah dilogaritmakan dengan sebelum dilogaritmakan adalah * Distribusinya, dengan menggunakan logaritma, kita mendekatkan distribusinya kepada distribusi normal. * Outliers-nya, dengan menggunakan logaritma, kita mengurangi banyak outliers yang ada pada kumpulan data tersebut.

Kita tertarik dalam mengestimasi koefisien dan intercept mereka, dalam nilai harga yang sesuai dan efek marginal dari peningkatan luas tanah berdasarkan harga.

mod3 <- lm(log(price)~sqft, data=br)
mod3$coefficients
##  (Intercept)         sqft 
## 1.083860e+01 4.112689e-04

Koefisiennya adalah \(b_1=10.84\) dan \(b_2=0.00041\), Menunjukan bahwa peningkatan luas tanah (sqft) dari apartemen setiap 1 unit (1 sqft) meningkatkan harga apartementnta sebesar 0.041 persen. Jadi, untuk rumah yang berharga $100,000, setiap peningkatan 100 sqft akan meningkatkan harga sekitar 100∗0.041 persen, yang bernilai $4112.7. Secara general, efek marginal dari peningkatan \(X\) pada \(Y\) dalam persamaan 5 adalah

\[ \begin{align} \tag{6} {dy\over dx}=\beta_2 * Y \end{align} \]

Dan untuk elastisitasnya adalah

\[ \begin{align} \tag{7} e= {dy\over dx} {X\over Y}=\beta_2 * X \end{align} \]

Beberapa perintah dibawah, menunjukan bagaimana cara menggambar nilai kurva dari model log-linear dan bagaimana cara menghitung efek marginal dan elastisitas untuk nilai tengah harga pada dataset. Nilai yang tepat dapat dihitung menggunakan formula berikut:

\[ \begin{align} \tag{8} \hat{Y}= e^{b_1+b_2*X} \end{align} \]

ordat <- br[order(br$sqft), ] #menurutkan dataset
mod4 <- lm(log(price)~sqft, data=ordat)
plot(br$sqft, br$price, col="grey")
lines(exp(fitted(mod4))~ordat$sqft, 
      col="blue", main="Log-linear Model")

pricex<- median(br$price)
sqftx <- (log(pricex)-coef(mod4)[[1]])/coef(mod4)[[2]]
(DyDx <- pricex*coef(mod4)[[2]])
## [1] 53.46495
(elasticity <- sqftx*coef(mod4)[[2]])
## [1] 0.9366934

R Memungkinkan kuta untuk menghitung kuantitas yang sama untuk beberapa (luas tanah, harga) pasangan dalam waktu yang sama, sepetri yang dilihatkan berikut:

b1 <- coef(mod4)[[1]]
b2 <- coef(mod4)[[2]]
#mengambul beberapa nilai untuk luas tanah(sqft):
sqftx <- c(2000, 3000, 4000) 
#mengestikasi harga untuk luas tanah dan menambahkan lagi:
pricex <- c(100000, exp(b1+b2*sqftx)) 
#menghitung ulang luas tanah (sqft) untuk seluruh harga:
sqftx <- (log(pricex)-b1)/b2 
#menghitung kan menununjukan elastisitasnya:
(elasticities <- b2*sqftx) 
## [1] 0.6743291 0.8225377 1.2338066 1.6450754

5 Variabel Indokator

Indikator atau variabel biner menunjukan wujud dari beberapa atribut dari unit yang di observasi, seperti jenis kelamin maupun suku jika unit yang di observasi bersifat individual atau lokasi jika unit yang di observasi adalah rumah. Dalam dataset utown, variabel utown adalah 1 jika rumah dekat dengan kampus dan 0 untuk sebaliknya. Ini adalah contoh model regresi linear yang melibatkan variabel utown:

\[ \begin{align} \tag{9} price_i= \beta_1+\beta_2*utwon_i \end{align} \]

Koefisien dari variabel dalam model linear sederhana sama dengan perbedaan rata-rata harga dari kedua kategori; koefisien intercept dari model dalam persamaan 9 sama dengan rata-rata harga dari rumah yang tidak dekat dengan kampus. Mari pertama-tama kita hitung rata-rata dari harga untuk setiap kategori, dimana dinyatakan dengan kode price0bar dan price1bar:

data(utown)
price0bar <- mean(utown$price[which(utown$utown==0)])
price1bar <- mean(utown$price[which(utown$utown==1)])

Hasil yang didapatkan adalah: \(\bar{price}=\) 277.2416012 dekat dengan kampus, dan \(\bar{price}=\) 215.7324948 untuk yang tidak dekat. Sekarang kita menunjukan bahwa hasil yang sama dihasilkan dari koefisien model regresi dalam Persamaan 9 :

mod4 <- lm(price~utown, data=utown)
b1 <- coef(mod4)[[1]] 
b2 <- coef(mod4)[[2]]
b1
## [1] 215.7325
b2
## [1] 61.50911

Hasil yang didapatkan adalah: \(\bar{price}=b_1=\) 215.7324948 untuk yang rumahnya jauh dari kampus, dan \(\bar{price}=b_1+b_2=\) 277.2416012 yang rumahnya dekat dari kampus.

6 Simulasi Monte Carlo

Simulasi Monte Carlo membuat nilai acak untuk variabel tetap ketika koefisien regresi dan distribusi random diberikan. Contoh berikut dilakukan untuk menentukan distribusi dari variabel bebas pada model pengeluaran makanan dalam persamaan 2,

N <- 40
x1 <- 10
x2 <- 20
b1 <- 100
b2 <- 10
mu <- 0
sig2e <- 2500
sde <- sqrt(sig2e)
yhat1 <- b1+b2*x1
yhat2 <- b1+b2*x2
curve(dnorm(x, mean=yhat1, sd=sde), 0, 500, col="blue")
curve(dnorm(x, yhat2, sde), 0,500, add=TRUE, col="red")
abline(v=yhat1, col="blue", lty=2)
abline(v=yhat2, col="red", lty=2)
legend("topright", legend=c("f(y|x=10)", 
                            "f(y|x=20)"), lty=1,
       col=c("blue", "red"))

Kemudian, kita menghitung varians dari \(b_2\) dan menunjukan plot dari density funtion berikut.

\[ \begin{align} \tag{10} var(b_2)={\sigma^2 \over \sum(x_i−\bar{x})} \end{align} \]

x <- c(rep(x1, N/2), rep(x2,N/2))
xbar <- mean(x)
sumx2 <- sum((x-xbar)^2)
varb2 <- sig2e/sumx2
sdb2 <- sqrt(varb2)
leftlim <- b2-3*sdb2
rightlim <- b2+3*sdb2
curve(dnorm(x, mean=b2, sd=sdb2), leftlim, rightlim)
abline(v=b2, lty=2)

Sekarang, dengan nilai \(b_1\), \(b_2\), dan error standar deviasi yang sama, kita dapat membuat set nilai untuk \(Y\) , regresi \(Y\) terhadap \(X\), dan menghitung nilai estimasi untuk koefisien \(b_2\)dan standar errornya.

set.seed(12345)
y <- b1+b2*x+rnorm(N, mean=0, sd=sde)
mod6 <- lm(y~x)
b1hat <- coef(mod6)[[1]]
b2hat <- coef(mod6)[[2]]
mod6summary <- summary(mod6) #summarynya terdapat nilai standard errors
seb2hat <- coef(mod6summary)[2,2]

Hasilnya adalah \(b_2=11.64\) dan \(s_e(b_2)=1.64\). Kelebihan dari Simulasi Monte Carlo adalah, kemungkinan dari pengulangan estimasi tingkat regresi untuk angka yang besar akan secara otomatis membuat sample. Jadi, kida dapat mendapatkan angka yang besar untuk parameter, katakanlah \(b_2\), dan kemudian memperkitakan karakteristik samplingnya. Sebagai contoh, jika rata-rata dari nilai tersebut mendekati nilai asumsi yaitu \(b_2= 10\), kita dapat menyimpulkan bahwa estimator (metode dari mengestimasi parameter) adalah tidak bias.

Sekarang kita akan menggunakan nilai dari \(X\) dalam dataset food, dan menghasilkan \(Y\) menggunakan model linear dengan \(b_1=100\) dan \(b_2= 10.\)

data("food")
N <- 40
sde <- 50
x <- food$income
nrsim <- 1000
b1 <- 100
b2 <- 10
vb2 <- numeric(nrsim) #menyimpan estimasi nilai b2
for (i in 1:nrsim){
  set.seed(12345+10*i)
  y <- b1+b2*x+rnorm(N, mean=0, sd=sde)
  mod7 <- lm(y~x)
  vb2[i] <- coef(mod7)[[2]]
}
mb2 <- mean(vb2)
seb2 <- sd(vb2)

Rata-rata dan standar deviasi dari estimasi 40 nilai dari \(b_2\) adalah, berturut-turut, 9.974985 dan 1.152632. Gambar dibawah menunjukan distribusi data yang telah disimulasikan dari \(b_2\) dengan yang dijelaskan secara teori.

plot(density(vb2))
curve(dnorm(x, mb2, seb2), col="red", add=TRUE)
legend("topright", legend=c("true", "simulated"), 
       lty=1, col=c("red", "black"))

hist(vb2, prob=TRUE, ylim=c(0,.4))
curve(dnorm(x, mean=mb2, sd=seb2), col="red", add=TRUE)

rm(list=ls()) # Peringatan: ini digunakan untuk membersihkan Environment