Permulaan/Pengertian

Regresi Linier Sederhana ini digunakan dalam memprediksi hasil kuantitatif y dengan berdasarkan satu variabel prediktor tunggal x, dimana bertujuan dalam membangun model atau rumus matematika yang mendefinisikan y sebagai fungsi dari variabel x, lalu setelah itu model yang sudah di signifikan secara statistik dapat digunakan dalam memprediksi hasil pada masa mendatang berdasarkan nilai x baru.

Dalam asumsi pemeriksan regresi, akan diperiksa sisa distribusi dimana regresi linier membuat beberapa asumsi tentang data, yaitu:

• Data Linieritas dimana merupakan hubungan antara prediktor (x) dan hasil linier (y)

• Residu Normalitas dimana sebagian besar residu terdistribusi normal.

• Homogenitas Varians Residual dimana residu memiliki varians yang konstant atau homoskedastisitas.

• Independensi dimana independensi ini merupakan suatu istilah kesalahan residual.

Diperlukan pemeriksaan apakah asumsi ini benar atau tidak. Dimana potensi masalah meliputi:

• Hasil Nonlinearitas yang merupakan hubungan prediktor

• Heteroskedastisitas yang merupakan sebuah varians istilah kesalahan yang tidak konstan

• Terdapat nilai yang berpengaruh dimana nilai-nilai berpengaruh dalam data yang dapat berupa:

  • Pencilan merupakan nilai ekstrim dalam variabel hasil (y)
  • Poin Leverage Tinggi merupakan nilai ekstrim dalam variabel prediktor (x).

Dimana semua asumsi dan masalah potensial ini diperiksa dengan membuat beberapa plot diagnostik yang memvisualisasikan sisa kesalahan.

Prediksi VS Peramalan

Mengenai asumsi bahwa nilai ekspektasi dari error term dalam Persamaan 1 adalah nol, Persamaan 2 memberikan y^i, nilai prediksi dari ekspektasi yi yang diberikan xi, di mana b1 dan b2 merupakan sebuah estimasi atau kuadrat terkecil dari parameter regresi β1 dan β2.

Nilai prediksi y^i merupakan variabel acak dikarenakan bertumpu pada sampel yang dapat menghitung interval kepercayaan dan menguji hipotesis tentangnya. Lalu juga dapat menentukan distribusi dan variansnya, dimana prediksi memiliki distribusi normal menjadi kombinasi linier dari dua variabel acak yang terdistribusi normal b1 dan b2, dan variansnya diberikan oleh Persamaan 3.

Ini diperlukan penggunaan estimasi varians agar dapat menggunakan distribusi t dan bukan distribusi normal. Persamaan 3 berlaku dalam setiap x yang diberikan, atau bisa disebut x0, tidak hanya untuk x - s dalam kumpulan data.

dimana kan dapat direduksi menjadi

Varians pada Persamaan 3 tidak sama dengan yang ada pada Persamaan 4.Pertama adalah varians dari perkiraan ekspektasi y, sedangkan yang terakhir merupakan varians dari kejadian tertentu y. Untuk itu persamaan 4 dapat disebut dengan varians dari kesalahan perkiraan. Varians kesalahan ramalan ini lebih besar dari varians dari prediksi E(y|x).

Dalam menentukan kesalahan standar untuk persamaan makanan untuk pendapatan rumah tangga (pendapatan) $2000 seminggu, yaitu pada persamaan x = x0 = 20 dengan menggunakan Persamaan 4 dengan mengambil var(b2) dan σ^, kesalahan standar regresi dari keluaran regresi.

library('remotes')
#install.packages("PoEdata")
library('PoEdata')
data("food")
alpha<-0.05
x<-20
xbar<-mean(food$income)
linmod<-lm(food_exp~income, data=food)
b1<-coef(linmod)[[1]]
b2<-coef(linmod)[[2]]
yhatx<-b1+b2*x
slinmod<-summary(linmod)
df<-df.residual(linmod)
tcr<-qt(1-alpha/2, df)
N<-nobs(linmod)
N<-NROW(food)
varb2<-vcov(linmod)[1,2]
sighat2<-slinmod$sigma^2
varf<-sighat2+sighat2/N+(x-xbar)^2*varb2
sef<-sqrt(varf)
lb<-yhatx-tcr*sef
ub<-yhatx+tcr*sef
print(c(lb,ub))

## [1] 104.2899 470.9279

Maka diperoleh hasil interval kepercayaan untuk ramalan [104.2899, 470.9279], dimana mengartikan bahwa lebih besar dari interval kepercayaan dari perkiraan nilai y seperti yang sudah diharapkan berdasarkan Persamaan nomor 4.Dari hasil tersebut dapat divisualisasikan dengan menggambar batas atas dan bawah bersama-sama dengan garis regresi. Gambar berikut menunjukkan pita interval kepercayaan tentang garis regresi.

sef<-sqrt(sighat2+sighat2/N+(food$income-xbar)^2*varb2)
yhatv<-fitted.values(linmod)
lbv<-yhatv-tcr*sef
ubv<-yhatv+tcr*sef
xincome<-food$income
dplot<-data.frame(xincome, yhatv, lbv, ubv)
dplotord<-dplot[order(xincome),]
xmax<-max(dplotord$xincome)
xmin<-min(dplotord$xincome)
ymax<-max(dplotord$ubv)
ymin<-min(dplotord$lbv)

plot(dplotord$xincome, dplotord$yhatv,
     xlim=c(xmin,xmax),
     ylim=c(ymin,ymax),
     xlab="Income", 
     ylab="Food Expenditure",
     type="1")
lines(dplotord$ubv~dplotord$xincome,lty=2)
lines(dplotord$lbv~dplotord$xincome,lty=2)

Cara lain dalam menemukan taksiran titik dan interval untuk prediksi E(y|x) dan perkiraan y adalah dengan menggunakan fungsi predict() pada R yang berfungsi dalam mensyaratkan nilai-nilai independen variabel di mana prediksi (atau ramalan) dimaksudkan memiliki struktur kerangka data. Contoh berikutnya menunjukkan pada titik paralel dan perkiraan interval pengeluaran makanan yang diprediksi dan diperkirakan untuk pendapatan adalah $2000. Seperti contoh diatas telah ditunjukkan bahwa perkiraan titik sama untuk prediksi dan ramalan hanya tetapi perkiraan interval sangat berbeda.

incomex=data/frame(income=20)
predict(linmod, newdata=incomex, interval="confidence", level=0.95)

predict(linmod, newdata=incomex, interval="prediction", level=0.95)

Menggunakan fungsi predict() untuk mereplikasi pita estimasi interval, dan titik-titik dalam kumpulan data.

xmax<-max(food$income)
xmin<-min(food$income)
income<-seq(from=xmin, to=xmax)
ypredict<-predict(linmod, newdata=data.frame(income), interval="confidence")
yforecast<-predict(linmod, newdata=data.frame(income), interval="prediction")
matplot(income, cbind(ypredict[,1], ypredict[,2], ypredict[,3], yforecast[,2], yforecast[,3]),
        type="1",lty=c(1,2,2,3,3),
        col=c("black","red","red","blue","blue"),
        ylab="Food Expenditure", xlab="Income")
points(food$income, food$food_exp)
legend("topleft",
       legend=c("E[y|x","lwr_pred","upr_pred","lwr_forcst","upr_forcst"),
       lty=c(1,2,2,3,3),
       col=c("black","red","red","blue","blue"))

Gambar tersebut menampilkan pita yang diprediksi dan yang diperkirakan pada grafik yang sama, untuk menunjukkan bahwa mereka memiliki perkiraan titik yang sama (garis hitam, solid) dan bahwa pita yang diperkirakan jauh lebih besar daripada yang diprediksi. Dengan kata lain, Anda mungkin berpikir tentang perbedaan antara dua jenis interval yang kami sebut prediksi dan ramalan sebagai berikut:

• Interval prediksi tidak seharusnya mencakup 95 persen poin, tetapi untuk memasukkan garis regresi E(y|x), dengan probabilitas 95 persen.

• Interval yang diperkirakan Pada sisi lain, harus mencakup setiap titik sebenarnya dengan probabilitas 95 persen.

Koefisien Determinasi

Analisis data betujuan dalam menjelaskan variabilitas dalam variabel dependen yaitu y. Untuk Membandingkan model, hal yang harus dilakukan adalahMempertimbangkan situasi dasar dimana yang dimiliki informasi tentang variabel dependen saja. Diperlukan pemodelan mean dari y, y^, untuk menjelaskan variabilitas, dimana dianggap sebagai model dasar atau baseline. Namun ketika memiliki informasi lain yang terkait dengan variabel dependen, seperti variabel independen kontinu, x maka dapat mencoba mengurangi variabilitas variabel dependen dengan menambahkan x ke model lalu memodelkan garis regresi linier yang dibentuk oleh nilai prediksi oleh x, y^. Untuk itu diharapkan dapat mengurangi variabilitas y sebagai hasil dari penambahan x ke dalam model. Oleh karena itu, perlu adanya evaluasi terhadap garis regresi, seperti seberapa jauh lebih baik daripada model dasar yang hanya menggunakan rata-rata variabel terikat.

Variasi total y terhadap mean sampelnya, Sum of Squares of Total (SST), dapat dijelaskan menjadi variasi terhadap garis regresi, Sum of Squares of Errors (SSE), dan variasi garis regresi terhadap mean dari y, Sum of Squares due to Segression (SSR ), seperti yang ditunjukkan oleh Persamaan 5.

Jika membutuhkan jumlah kesalahan kuadrat, SSE, atau jumlah kuadrat karena regresi, SSR, dapat menggunakan fungsi anova yang strukturnya ditunjukkan pada Tabel 1 berikut.

library(stats)
#install.packages("hms")
library(kableExtra)
library(DT)
anov<-anova(linmod)
SSE<-anov[2,2]
SSR<-anov[1,2]
SST<-anov[1,2]+anov[2,2]
dfr<-data.frame(anov)
datatable(dfr,
          caption=htmltools::tags$caption(
            style='caption-side: bottom; text-align:center;',
            htmltools::em('Table 1:Output generated by the anova function,')),
          options=list(dom='t'))

Hasil ANOVA pada Tabel 1 berisi informasi berguna lainnya: jumlah derajat kebebasan, anov[2,1] dan varians yang diperkirakan sigma kudrat = anov[2,3]. Hasil ini juga menjelaskan, koefisien determinasi, R2, didefinisikan sebagai proporsi varians dalam y yang dijelaskan oleh regresi, SSR, dalam variasi total y, SST. Membagi kedua ruas Persamaan 5 dengan SST dan menyusun suku-sukunya memberikan rumus untuk menghitung R2, seperti yang ditunjukkan pada Persamaan 6.

R2 bernilai antara 0 dan 1, dengan nilai yang lebih tinggi menunjukkan kecocokan yang lebih dekat dari garis regresi ke data (nilai R-kuadrat yang besar menunjukkan kecocokan yang baik). Dalam R, nilai R2 dapat diambil dari ringkasan model regresi dengan nama adj.r.squared misalnya, dalam contoh makanan diatas bahwa R^2=0,3688. Dimana hal itu juga dicetak sebagai bagian dari ringkasan model regresi, seperti yang ditunjukkan padaurutan kode berikut.

slinmod

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

(rsq<-slinmod$adj.r.squared)

## [1] 0.3688181

Dalam model regresi sederhana ini jumlah kuadrat karena regresi hanya mencakup pendapatan variabel. Dalam contoh makanan kita, R2=0,3688. yang menunjukkan model tersebut tidak cocok dengan data yang diamati. Untuk memastikan apakah ada hubungan linier yang signifikan secara statistik antara variabel dependen dan independen maka akan dilakukan dengan mempertimbangkan uji T: uji signifikansi kemiringan. Lalu, dalam model regresi berganda yaitu model dengan lebih dari satu variabel bebas, jumlah kuadrat akibat regresi sama dengan jumlah kuadrat semua variabel bebas, akan menggunakan uji F Keseluruhan: evaluasi model global.

Evaluasi Kemiringan Individual

Pada regresi linier sederhana dengan hanya satu x, hasil uji-F global dan Uji-T (signifikansi kemiringan) memiliki kesimpulan yang sama. Uji signifikansi formal kemiringan adalah sebagai berikut:

Statistik uji: \(t=b1/SE_(b1)\), di mana SE_(b1) merupakan kesalahan standar b1. Tes ini menanyakan apakah terdapat hubungan linier yang signifikan secara statistik antara variabel dependen dan independen. Jika tidak, maka model regresi tidak baik karena asumsi dasarnya adalah hubungan linier antara dua variabel. Jika tidak ada hubungan linier, maka kemiringan mendekati nol dan scatterplot y dan x diharapkan menunjukkan sebaran acak, yang berarti tidak ada korelasi.

Dengan melakukannya, akan mendapatkan nilai p yang sesuai dengan nilai t yang dihitung dan merujuk ke situs web’nilai p dari score calculator atau score calculator F` . Dalam hal ini,diberikan ‘skor T’ sebagai 4,877, dan ‘DF’ sebagai 38, yang merupakan df dari MSE, n-2, dan akhirnya mendapatkan nilai p sebagai 1,946e-05, yang merupakan nilai p yang sama dari F-test di atas.

Jika ingin menemukan nilai-p pada R, dapat dijalankan kode berikut.

anova(linmod)$'Pr(>F)`[1]

Dari hasil tersebut diketahui bahwa p-value = 1.945862e-05 lebih kecil dari taraf signifikansi =0.05. Oleh karena itu, menolak hipotesis nol bahwa kemiringannya nol lalu dapat disimpulkan bahwa Terdapat hubungan linier yang signifikan antara y dan x.

Nilai yang Pas dan Residual

Model regresi linier membuat asumsi bahwa hubungan antara prediktor (x) dan variabel hasil adalah linier. Ini mungkin tidak benar. Hubungannya bisa polinomial atau logaritmik. Selain itu, data mungkin berisi beberapa pengamatan yang berpengaruh, seperti outlier (atau nilai ekstrim), yang dapat mempengaruhi hasil regresi. Oleh karena itu, diharuskan mendiagnosa model regresi yang dibangun dalam mendeteksi potensi masalah dan memeriksa apakah asumsi yang dibuat oleh model regresi linier terpenuhi atau tidak.

Dengan begitu dilakukan pemeriksaan distribusi kesalahan residual, yang dapat menginformasikannya lebih banyak tentang data yang akan dianalisis. Pada R dapat dengan mudah menambah data yang betujuan untuk menambahkan nilai dan residu yang sesuai dengan menggunakan fungsi augment() [broom package]. atau biasa disebut model keluaran model.diag.metrics karena berisi beberapa metrik yang berguna dalam mendiagnostik regresi.

#install.packages("broom")
library(broom)
library(DT)
model.diag.metrics<-augment(linmod)
datatable(model.diag.metrics,
          caption=htmltools::tags$caption(
            style='caption-side: bottom; text-align:center;',
            htmltools::em('Table 2:The distriburion of Residuals errors.')),
          extensions='FixedColumns',
          options=list(scrollX=TRUE, fixedColumns=TRUE))

Langkah selanjutnya adalah memvisualisasikan kesalahan residual (dalam warna merah) antara nilai yang diamati dan garis regresi yang dipasang. Setiap segmen merah vertikal mewakili kesalahan residual antara nilai penjualan yang diamati dan nilai prediksi yang pas dan sesuai.

library(ggplot2)
ggplot(model.diag.metrics, aes(income, food_exp))+
  geom_point()+
  stat_smooth(method=lm, se=FALSE)+
  geom_segment(aes(xend=income, yend=.fitted),color="red",size=0.3)

Diagnosa Regresi

Pada Statistik, Diagnostik Regresi merupakan salah satu dari bebeerapa prosedur yang tersedia dalam menganalisis regresi yang berusaha menilai validitas model dengan dalam cara yang berbeda. Plot diagnostik menunjukkan residu dalam empat cara berbeda:

• Residuals vs Fitted Digunakan dalam memeriksa asumsi hubungan linier. Garis horizontal, tanpa pola yang berbeda merupakan indikasi untuk hubungan linier.

• Normal Q-Q Digunakan dalam memeriksa apakah residual berdistribusi normal. Namun, ada baiknya jika titik residu mengikuti garis putus-putus lurus.

• Scale-Location Digunakan dalam memeriksa homogenitas varians dari residual (homoscedasticity), dimana garis horizontal dengan titik penyebaran yang sama merupakan indikasi yang baik dari homoskedastisitas.

• Residuals vs Leverage Digunakan dalam mengidentifikasi kasus-kasus yang berpengaruh, yaitu nilai-nilai ekstrim yang mungkin mempengaruhi hasil regresi ketika dimasukkan atau dikeluarkan dari analisis.

Plot diagnostik dapat dibuat menggunakan fungsi dasar R plot() atau fungsi autoplot() [paket ggfortify], yang membuat grafik berbasis ggplot2:

library(dplyr)
#install.packages("ggfortify")
library(ggfortify)
par(mfrow=c(2,2))
autoplot(linmod)

Di bagian berikut akan menjelaskan secara rinci cara menggunakan grafik dan metrik ini untuk memeriksa asumsi regresi dan untuk mendiagnosis masalah potensial di masa depan.

set.seed(12345)
x<-runif(1000,-2.5,2.5)
e<-rnorm(1000,0,4)
y<-15-4*x^2+e
quadmod<-lm(y~x)
ehat<-resid(quadmod)
ymi<-min(ehat)
yma<-max(ehat)
plot(x, ehat, ylim=c(ymi,yma),xlab="x",ylab="residuals",col="grey")

#install.packages("quantmod")
library(quantmod)
#install.packages('tseries')
library(tseries)
ehat<-resid(linmod)
ebar<-mean(ehat)
sde<-sd(ehat)
hist(ehat, col="grey",freq=FALSE, main="",ylab="density",xlab="ehat")
curve(dnorm(x,ebar,sde),col=2,add=TRUE,ylab="density",xlab="ehat")

jarque.bera.test(ehat)

## 
##  Jarque Bera Test
## 
## data:  ehat
## X-squared = 0.06334, df = 2, p-value = 0.9688

linlogmod<-lm(log(food_exp)~income, data=food)
plot(linlogmod,3)

high_lev<-data.frame(
  x=c(food$income,500,700),
  y=c(food$food_exp,8000,10000)
)
high_lev_mod<-lm(y~x,high_lev)
par(mfrow=c(1,2))
plot(high_lev_mod,4)
plot(high_lev_mod,5)

Model Kuadratik

Model pada bagian sebelumnya terlihat bagus, tetapi dapat melihat bahwa plot memiliki kelengkungan yang tidak dijelaskan dengan baik oleh model linier. Sekarang cocok dengan model pendapatan kuadrat sehingga akan membuat variabel yang disebut income^2 yang merupakan kuadrat dari pendapatan variabel.

#install.packages("remotes")
remotes::install_github("ccolonescu/PoEdata")

library(PoEdata)
data("food")
food$income2<-food$income^2
quadmod1<-lm(food_exp~income2, data=food)
summary(quadmod1)

## 
## Call:
## lm(formula = food_exp ~ income2, data = food)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -219.142  -56.559   -6.527   72.755  198.143 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 166.45749   29.22934   5.695 1.49e-06 ***
## income2       0.27232    0.05907   4.610 4.46e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 91.41 on 38 degrees of freedom
## Multiple R-squared:  0.3587, Adjusted R-squared:  0.3418 
## F-statistic: 21.25 on 1 and 38 DF,  p-value: 4.455e-05

Sintaks yang terlibat dalam pemasangan model linier dengan dua atau lebih prediktor. Kami menyertakan setiap prediktor dan memberi tanda plus di antara mereka.

Hasil dari model kuadratik di atas ternyata tidak sesuai dengan yang diharapkan. langkah selanjutnya adalah melihat contoh yang baik dalam model kuadrat dengan memahami hubungan antara jumlah jam kerja dan kebagiaan. Dimana kita memiliki data berikut tentang jumlah jam kerja per minggu dan tingkat kebahagiaan yang dilaporkan (pada skala 0-100) untuk 11 orang yang berbeda:

data<-data.frame(hours=c(6,9,12,14,30,35,40,47,51,55,60),
                 happiness=c(14,28,50,70,89,94,90,75,59,44,27))
data$hours2<-data$hours^2
quadmod2<-lm(happiness~hours+hours^2, data=data)
summary(quadmod2)

## 
## Call:
## lm(formula = happiness ~ hours + hours^2, data = data)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -39.34 -21.99  -2.03  23.50  35.11 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  48.4531    17.3288   2.796   0.0208 *
## hours         0.2981     0.4599   0.648   0.5331  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 28.72 on 9 degrees of freedom
## Multiple R-squared:  0.0446, Adjusted R-squared:  -0.06156 
## F-statistic: 0.4201 on 1 and 9 DF,  p-value: 0.5331

Model Polinomial

Contoh penerapan model polinomial pada data food.

food$income2<-food$income^2
polymod1<-lm(food_exp~income+I(income2), data=food)
summary(polymod1)

## 
## Call:
## lm(formula = food_exp ~ income + I(income2), data = food)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -223.037  -49.992   -5.946   67.605  212.458 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) 81.982928  73.089590   1.122    0.269
## income      10.405746   8.263941   1.259    0.216
## I(income2)  -0.005607   0.228375  -0.025    0.981
## 
## Residual standard error: 90.72 on 37 degrees of freedom
## Multiple R-squared:  0.385,  Adjusted R-squared:  0.3518 
## F-statistic: 11.58 on 2 and 37 DF,  p-value: 0.0001242

Dari ringkasan polinomial tersebut dapat dilihat bahwa model tidak sesuai dengan yang diharapkan. Oleh karena itu, model iini tidak bisa diterapkan dalam kumpulan data makanan..

Selanjutnya, terdapat model regresi lain yang mungkin menyertakan istilah kuadrat atau kubik untuk lebih menggambarkan sifat data. Berikut ini adalah contoh model kuadrat dan kubik untuk dataset wa_wheat, yang memberikan hasil gandum tahunan dalam ton per hektar di Greenough Shire di Australia Barat selama periode 48 tahun.

library(PoEdata)
data("wa_wheat")
wa_wheat$time<-wa_wheat$time^2
polymod2<-lm(greenough~time+I(time^2),data=wa_wheat)
summary(polymod2)

## 
## Call:
## lm(formula = greenough ~ time + I(time^2), data = wa_wheat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.47116 -0.15002  0.00728  0.16748  0.31939 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 8.719e-01  5.069e-02  17.200   <2e-16 ***
## time        1.572e-04  1.339e-04   1.173   0.2468    
## I(time^2)   1.402e-07  6.347e-08   2.209   0.0323 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1848 on 45 degrees of freedom
## Multiple R-squared:  0.7551, Adjusted R-squared:  0.7442 
## F-statistic: 69.39 on 2 and 45 DF,  p-value: 1.782e-14

Model Linear-Log

Bentuk fungsional non-linier dari model regresi berguna untuk hubungan antara dua variabel yang lebih kompleks daripada yang linier. Solusi yang mungkin untuk memecahkan masalah fungsional ini adalah dengan menggunakan bentuk non-linier berdasarkan model matematika, penalaran, atau hanya memeriksa plot pencar dari data. Dalam model pengeluaran makanan, misalnya, masuk akal untuk meyakini bahwa jumlah yang dibelanjakan untuk makanan meningkat lebih rendah daripada pendapatan yang lebih tinggi. Fungsi logaritmik cocok dengan profil ini dan ternyata, relatif mudah untuk diraih, yang sangat populer dalam model ekonometrik.

Bentuk umum dari model ekonometrika log linier diberikan pada Persamaan 8.

Efek marjinal dari perubahan x pada y adalah kemiringan kurva regresi dan diberikan oleh Persamaan 11, tidak seperti dalam bentuk linier,dimana itu bergantung pada x dan oleh karena itu, hanya berlaku untuk perubahan kecil di x.

Terkait dengan model log linier, dapat mengukur yang lain dimana menarik dalam ekonomi adalah semi-elastisitas y terhadap x, yang diberikan oleh Persamaan 12 . Semi-elastisitas menunjukkan bahwa perubahan x sebesar 1% mengubah y sebesar 2/100% unit y. Karena semi-elastisitas juga berubah ketika x berubah, itu seharusnya hanya ditentukan untuk perubahan kecil pada x.

Besaran lain yang mungkin menarik adalah elastisitas y terhadap x, yang diberikan oleh Persamaan 13 dan menunjukkan bahwa satu persen peningkatan x menghasilkan (β2/y) persen perubahan y.

Langkah selanjutnya adalah memperkirakan model log linier dalam dataset makanan, menggambar kurva regresi, dan menghitung efek marjinal untuk beberapa nilai tertentu dari variabel dependen.

#install.packages("xtable")
library(xtable)
linlogmod<-lm(food_exp~log(income), data=food)
tbl<-data.frame(xtable(linlogmod))
datatable(tbl,
          caption=htmltools::tags$caption(
            style='caption-side: bottom; text-align: center;',
            htmltools::em('Table 3: Linear-Log Model Output for the food example.')),
          options=list(dom='t')
          )

b1<-coef(linlogmod)[[1]]
b2<-coef(linlogmod)[[2]]
pred_linlogmod<-predict(linlogmod,newdata=data.frame(income),
                        interval="confidence")
plot(food$income, food$food_exp, xlab="Income",ylab="Food Expendicture")
lines(pred_linlogmod[,1]~income,lty=1,col="black")
lines(pred_linlogmod[,2]~income,lty=2,col="red")
lines(pred_linlogmod[,3]~income,lty=2,col="red")

x<-10
y<-b1+b2*log(x)
DyDx<-b2/x
DyPDx<-b2/100
PDyPDx<-b2/y
print(c(DyDx, DyPDx, PDyPDx))

## [1] 1.02096430 0.10209643 0.09548453

Hasil dalam pendapatan sebesar \(1000\) adalah sebagai berikut: \(dy/dx=13.217\), yang menunjukkan bahwa peningkatan pendapatan sebesar \(100\) (yaitu, satu unit x )untuk meningkatkan pengeluaran sebesar \(13.217\) ; untuk peningkatan pendapatan sebesar 1%, yaitu peningkatan sebesar \(10\), pengeluaran meningkat sebesar \(1,322\) dan, akhirnya, untuk peningkatan 1% dalam pengeluaran pendapatan meningkat sebesar 0,638%.

Model Log-Linear

Mengubah variabel tergantung dengan fungsi log() berguna ketika variabel memiliki distribusi miring, yang umumnya terjadi dengan jumlah yang tidak boleh negatif. Transformasi log() sering membuat distribusi mendekati normal. Model log-linier umum diberikan dalam Persamaan 14.

Rumus yang selanjutnya dengan mudah diturunkan dari persamaan log-linier 14. Semi-elastisitas di sini memiliki interpretasi yang berbeda dari yang ada dalam model log-linier: di sini, peningkatan x satu unit (dari x) menghasilkan perubahan 100b2 persen dalam y. Untuk perubahan kecil pada x, jumlah 100b2 dalam model log-linier juga dapat diinterpretasikan sebagai laju pertumbuhan dalam y (sesuai dengan peningkatan unit dalam x). Misalnya, jika x adalah waktu, maka 100b2 adalah laju pertumbuhan dalam y per satuan waktu.

• Prediksi \(y^n=exp(b1+b2x)\) atau \(exp(b1+b2x+σ^22)\) dengan prediktor “alami” y^n untuk digunakan dalam sampel kecil dan prediktor “dikoreksi”, y^c, dalam sampel besar

• Slope / efek marjinal : \(dy/dx=b_2y\)

• Semi-elastisitas : \(%Δy=100b_2Δx\)

library(xtable)
library(kableExtra)
loglinmod1 <- lm(log(greenough)~time, data=wa_wheat)
sloglinmod1 <-summary(loglinmod1)
tbl <- data.frame (xtable(sloglinmod1))
datatable(tbl, caption= htmltools::tags$caption(
  style='caption-side: botton; text-align: center;',
  htmltools::em("Table 4: Log-linear model for the yield equation.")),
  options=list (dom= 't'))

Tabel 4 memberikan b_2=0,017844 , yang menunjukkan bahwa laju pertumbuhan produksi gandum telah meningkat dengan laju rata-rata sekitar 1,78 persen per tahun.

Persamaan log-linear upah memberikan contoh lain untuk menghitung tingkat pertumbuhan, tetapi kali ini variabel bebasnya bukan waktu , tetapi pendidikan . Prediksi dan kemiringan dihitung untuk educ=12 tahun.

data("cps4_small", package="PoEdata")
xeduc=12
loglinmod2 <- lm(log(wage)~educ, data=cps4_small)
sloglinmod2 <-summary(loglinmod2)
tabl <- data.frame(xtable(sloglinmod2))
datatable(tabl,
caption=htmltools::tags$caption( 
style='caption-side: bottom, text-align: center;',
htmltools::em('Table 5: Log-linear wage regression output.')), 
options=list(dom='t'))

b1 <- coef(sloglinmod2)[[1]]
b2 <- coef(sloglinmod2)[[2]]
sighat2 <- sloglinmod2$sigma^2
g<- 100*b2                     # growth rate
yhatn<-exp(b1+b2*xeduc)        # "natural" predictiction
yhatc <- exp(b1+b2*xeduc+sighat2/2)# corrected prediction
DyDx <- b2*yhatn               #marginal effect
print(c(g,yhatn, yhatc, DyDx)) #print the result

## [1]  9.040825 14.795801 16.996428  1.337662

Berikut adalah hasil dari perhitungan tersebut: prediksi “alami” y^n=14.796. Prediksi terkoreksi, y^c=16.996 tingkat pertumbuhan g=9.041; dan efek marginal dy/dx=1.34. Laju pertumbuhan tersebut menunjukkan bahwa peningkatan pendidikan sebesar satu satuan (lihat uraian data menggunakan ?cps4_small) meningkatkan upah per jam sebesar 9.041 persen.

Gambar di bawah menyajikan garis regresi “alami” dan “dikoreksi” untuk persamaan upah, bersama dengan titik data aktual.

education=seq(0,22,2)
yn <- exp(b1+b2*education)
yc <- exp(b1+b2*education+sighat2/2)
plot(cps4_small$educ, cps4_small$wage, 
     xlab="education", ylab="wage",col="grey") 
lines (yn~education, lty=2, col="black")
lines(yc~education, lty=1, col="blue") 
legend("topleft", legend=c("yc","yn"),lty=c(1,2),col=c("blue","black"))

R^2 reguler tidak dapat digunakan dalamdua model regresi yang memiliki catatan variabel yang berbeda seperti model log-linier dan log-linier. Ketika itu membandingkan perbandingan seperti itu sangat diperlukan, seseorang dapat menggunakan R^2 umum, yaitu \(R^2g=[corr(y,y^]2\). Selanjutnya adalah menghitung R^2 umum untuk model kuadratik dan log linier

quadmod3 <- lm(wage~I(educ^2), data=cps4_small)
yhat3 <- predict (quadmod3) 
loglinmod3 <- lm(log(wage)~educ, data=cps4_small)
sloglinmod3 <- summary(loglinmod3)
b1 <- coef(sloglinmod3)[[1]] 
b2 <- coef(sloglinmod3)[[2]]
sighat2 <- sloglinmod3$sigma^2
yhat4 <- exp(b1+b2*cps4_small$educ+sighat2/2) 
rgquadmod3 <- cor(cps4_small$wage, yhat3)^2
rgloglinmod3 <- cor(cps4_small$wage, yhat4)^2
print(c(rgquadmod3, rgloglinmod3))

## [1] 0.1881762 0.1859307

Model kuadrat menghasilkan R^{2=0,188, dan model log-linier menghasilkan R}2=0,186. Dikarenakan yang pertama lebih tinggi, Lalu dapat disimpulkan bahwa model kuadrat lebih cocok dengan data daripada yang log-linear. (Namun, pengujian lain tentang bagaimana kedua model memenuhi asumsi regresi linier dapat mencapai kesimpulan yang berbeda; R^2 hanyalah salah satu kriteria pemilihan model.)

Untuk menentukan perkiraan interval perkiraan dalam model log-linier, pertama-tama dengan membangun interval dalam log menggunakan prediktor alami y^n, kemudian mengambil antilog dari batas interval. Peramalan sama seperti sebelumnya, diberikan dalam Persamaan 4. Perhitungan berikut menggunakan tingkat pendidikan sama dengan 12 dan = 0,05.

# The wage Log-Linear model 
#Prediction interval for educ= 12

alpha <- 0.05 
xeduc <- 12
xedbar <- mean(cps4_small$educ)
loglinmod4<- lm(log(wage)~educ, data=cps4_small)
b1 <- coef(loglinmod4)[[1]]
b2 <- coef(loglinmod4)[[2]]
df5 <- loglinmod4$df.residual
N <- nobs (loglinmod4) 
tcr <- qt(1-alpha/2, df=df5)
sloglinmod4 <- summary(loglinmod4)
varb2 <- vcov(loglinmod4) [2,2] 
sighat2 <- sloglinmod4$sigma^2
varf <- sighat2+sighat2/N+(xeduc-xedbar)^2*varb2
sef <- sqrt(varf) 
lnyhat <- b1+b2*xeduc
lowb <- exp(lnyhat-tcr*sef)
upb <- exp(lnyhat+tcr*sef) 
print(c(lowb, upb))

## [1]  5.260398 41.615811

Hasilnya adalah selang kepercayaan (5.26,41.62) . Selanjutnya, dengan pertimbangkan pita kepercayaan 95% untuk model upah log-linear di bawah ini:

#Drawing a confidence band for the Log-Linear
#wage equation

xmin <- min(cps4_small$educ) 
xmax <- max(cps4_small$educ)+2
education <- seq(xmin, xmax, 2)
lnyhat <- b1+b2*education
yhat <- exp(lnyhat) 
varf <- sighat2+sighat2/N+(education-xedbar)^2*varb2
sef <- sqrt(varf)
lowb <- exp(lnyhat-tcr*sef) 
upb <- exp(lnyhat+tcr*sef)
plot (cps4_small$educ, cps4_small$wage, col="grey", 
      xlab="education", ylab="wage", ylim=c(0,100))
lines (yhat~education, lty=1, col="black")
lines (lowb~education, lty=2, col="blue")
lines (upb~education, lty=2, col="blue") 
legend("topleft", legend=c("yhat", "lowb", "upb"),
       lty=c(1, 2, 2), col=c("black", "blue", "blue"))

Model Log-Log

Model log-log memiliki sifat yang diinginkan yaitu koefisien variabel bebas sama dengan elastisitas (konstan) y terhadap x. Oleh karena itu, model ini sering digunakan untuk mengestimasi persamaan penawaran dan permintaan. Bentuk standarnya diberikan dalam persamaan 15, di mana y,x, dan e adalah vektor N×1.

# Calculating Log-Log demand for chicken

data("newbroiler", package="PoEdata")
loglogmod <- lm(log(q)~log(p), data=newbroiler)
b1 <- coef(loglogmod)[[1]]
b2 <- coef(loglogmod)[[2]]
sloglogmod <- summary(loglogmod)
tbl <- data.frame(xtable (sloglogmod)) 
datatable (tbl, 
           caption=htmltools::tags$caption(
           style='caption-side: bottom; text-align: center;',
           htmltools::em ("Table 6: The log-log poultry regression equation.")),
           options=list(dom='t'))

Tabel 6 memberikan output regresi log-log. Koefisien pada p menunjukkan bahwa kenaikan harga sebesar 1% mengubah jumlah yang diminta sebesar 1.121% .

#Drawing the fitted values of the Log-Log equation

ngrid <-20 
xmin <- min(newbroiler$p)
xmax <-max(newbroiler$p)
step <- (xmax-xmin)/ngrid 
xp <- seq(xmin, xmax, step)
sighat2 <- sloglogmod$sigma^2
yhatc <- exp(b1+b2*log(newbroiler$p)+sighat2/2) 
yc <- exp(b1+b2*log(xp)+sighat2/2)
plot(newbroiler$p, newbroiler$q, ylim=c(10,60), xlab="price",ylab="quantity") 
lines(yc~xp, lty=1, col="black")

rgsq<-cor(newbroiler$q, yhatc)^2
rgsq

## [1] 0.8817758

R^2 umum dimana menggunakan nilai pas yang dikoreksi sama dengan 0,8818.

Referensi

• https://rpubs.com/dsciencelabs/econometrics4

• Bruce, Peter, and Andrew Bruce. 2017. Practical Statistics for Data Scientists. O’Reilly Media.

• James, Gareth, Daniela Witten, Trevor Hastie, and Robert Tibshirani. 2014. An Introduction to Statistical Learning: With Applications in R. Springer Publishing Company, Incorporated.