##R ile Çoklu Regresyon Analizi: Wooldridge Örnekleri ile
###İki bağımsız değişkenli model Bu bölümün Wooldridge kitabınızda bulunan ilk örnek maaş (wage) bağımlı değişkenin iki bağımsız değişkenle açıklamaya çalıştığımız model. Bu modeli tekrar hatırlayalım.
wage=β0+β1educ+β2exper+u
Burada exper iş hayatındaki deneyim yılı,educ eğitim yılını gösteren bağımsız değişkenler, ve u ise gözlemleyemediğimiz diğer faktörler olarak özetlenebilir .Burada basit regresyon modelinde de anlamaya çalıştığımız gibi asıl ilgilendiğimiz parametre β1 yani maaşı etkileyen diğer bütün faktörleri sabit tuttuğumuzda eğitimin gelir üzerinde ne kadar bir etkisi olduğu. Çoklu regresyon modelini kullanmamızın bir diğer nedeni de örneğimizden anlaşılacağı üzere gözlemlenemeyen (atlanan) değişken sapmasını ortadan kaldırabilmektir. Diğer bir örnek: Öğrenci başına yapılan harcamanın expend,liseöğrencilerinin ortalama test sonuçları avgscore üzerideki etkisini inceliyor.
avgscore=β0+β1expend+β2avginc+u
avginc ailenin ortalama geliri ve u diğer gözlemlenemeyen faktörleri temsil eder.Bulmaya çalıştığımız etki ,öğrenci başına yapılan expend ,lise öğrencilerinin ortalama test sonuçları avgscore üzerindeki etkisini gösteren b1’dir.Ancak ailenin ortalama geliri ,harcamalar üzerinde bir etkiye sahip olduğundan ,avginc değişkenini regresyona dahil etmemek ,en küçük kareler metoduyla bulduğumuz basit regresyon katsayısını b1 sapmalı hale getirecektir. İki bağımsız değişkenli modeli genel bir şekilde yazarsak :
y=β0+β1x1+β2x2+u
burada b0 kesen,b1 diğer faktörleri sabit tutarak y’deki değişimi x1’e göre ölçer b2 diğer faktörleri sabit tutarak y’deki değişimi x2’ye göre ölçer. İki açıklayıcı değişkenli regresyon modelinin kitabınızdaki bir diğer örneği, doğrusal regresyon modelinin fonksiyonel olarak yazılabildiğini göstermek için verilmiştir. Ailenin tüketimini, ailenin gelirinin ikinci dereceden (quadratic) bir fonksiyonu olarak gösterirsek:
cons=β0+β1inc+β2inc2+u
cons ailenin tüketimi,inc ailenin geliri ,u tüketimi etkileyen diğer faktörler olarak özetlenebilir. ##K BAGIMSIZ DEGIŞKENLI MODEL: Çoklu regresyon modeli olarak adlandırılabilir.
y=β0+β1x1+β2x2+β3x3+…+βkxk+u
u terimi hata terimi veya gözlenemeyen diğer faktörler olarak adlandırılabilir.
Çoklu regresyon modelinin temel varsayımını hatırlayalım.
E(u|x1,x2,x3,…,xk)=0
Gözleyemediğimiz diğer faktörlerin, bağımsız değişkenlerle ilişkisiz olması gerektiğini gösteren varsayımdır. ##ORNEK 3.1 KOLEJ GPS’SINI BELIRLEYICILERI
library(wooldridge)
data(gpa1)
library(rmarkdown)
paged_table(gpa1)
Örnek için kullanacağımız değişkenlerin tanımlarını şöyle yapabiliriz.
colGPA: MSU GPA, Michigan State Üniversitesi not ortalaması
hsGPA: high school GPA, lise not ortalaması
ACT: ‘achievement’ score, ACT, bir öğrencinin okulda ne öğrendiğini ölçen bir başarı testidir. SAT daha çok bir yetenek testidir, muhakeme ve sözel yetenekleri test eder. Üniversiteler SAT veya ACT sonuçları ile öğrenci kabul eder.
coklureg <- lm(colGPA~ hsGPA+ACT,data = gpa1)
summary(coklureg)
##
## Call:
## lm(formula = colGPA ~ hsGPA + ACT, data = gpa1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.85442 -0.24666 -0.02614 0.28127 0.85357
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.286328 0.340822 3.774 0.000238 ***
## hsGPA 0.453456 0.095813 4.733 5.42e-06 ***
## ACT 0.009426 0.010777 0.875 0.383297
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3403 on 138 degrees of freedom
## Multiple R-squared: 0.1764, Adjusted R-squared: 0.1645
## F-statistic: 14.78 on 2 and 138 DF, p-value: 1.526e-06
colGPA=1.29+0.453hsGPA+0.0094ACT+u
İlk olarak kesişim katsayısı, hsGPA ve ACT’nin her ikisi de sıfır olursa, öngörülen üniversite GPA’sini göstermektedir ve örneğimizde 1.29 olarak hesaplanmıştır. Üniversiteye giden hiç kimsenin lise not ortalaması sıfır veya başarı testi sıfır olmadığı için, bu denklemdeki kesişim kendi başına anlamlı değildir. Eğim katsayılarına bakınca ikisinin de bekleneceği gibi pozitif olduğunu görüyoruz. Eğer ACT puanını sabit tutarsak, hsGPA’deki 1 birimlik artış colGPA üzerinde 0.453 puanlık bir artış yaratacaktır. Başka bir deyişle, eğer ACT puanları aynı olan iki öğrenci varsa ve bir öğrencinin lise not ortalaması diğer öğrenciden bir puan yüksekse, bu öğrencinin üniversite not ortalamasının yarım puan daha fazla olmasını tahmin edebiliriz.
Üniversite not ortalamalarının 4 üzerinden, ACT sonuçlarının 36 puan üzerinden hesaplandığı düşünülürse, ACT’nin eğimi 0.0094’ün üniversite not ortalaması üzerinde fazla bir etkisi olmadığı düşünülebilir. Hata paylarını göz önüne aldığımızda ACT’nin sadece küçük değil aynı zamanda istatiksel olarak da anlamlı olmadığını göreceğiz.
Basit regresyon örneğimizi hatırlarsak:
basitreg <- lm(colGPA~ACT,data = gpa1)
summary(basitreg)
##
## Call:
## lm(formula = colGPA ~ ACT, data = gpa1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.85251 -0.25251 -0.04426 0.26400 0.89336
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.40298 0.26420 9.095 8.8e-16 ***
## ACT 0.02706 0.01086 2.491 0.0139 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3656 on 139 degrees of freedom
## Multiple R-squared: 0.04275, Adjusted R-squared: 0.03586
## F-statistic: 6.207 on 1 and 139 DF, p-value: 0.0139
ACT’nin katsayısı çoklu regresyona göre bir hayli büyük kalıyor (0.02706/0.009426=2.87 kat daha fazla). ACT değerinin çoklu regresyonda lise not ortalaması eklenince azaldığını gözlemliyoruz. Çünkü lise not ortalaması ve ACT arasında yüksek bir korelasyon var. Bu durum bize aslında basit regresyon kullanıldığında temel varsayımı ihlal ettiğimize dair güçlü bir ipucu veriyor.
Bu örnekte son olarak basit ve iki bağımsız değişkenli regresyon grafiklerini gösterelim. Basit regresyon örneğimiz için, dağılım grafiğine bir regresyon çizgisi eklememizin yeterli olduğunu hatırlayalım.
plot(gpa1$ACT,gpa1$colGPA)
abline(basitreg)
Bu grafik sadece iki boyutlu olduğundan, iki boyutlu bilgisayar ekranlarımızda göstermek daha kolaydır. Y eksenine üniversite not ortalaması, X eksenine ACT puanlarını koyuyoruz. Son olarak bulduğumuz eğim ve kesen katsayılarını grafiğimizde temsil eden doğrumuzu çiziyoruz. Ancak çok değişkenli regresyon grafiklerini, iki boyutlu bilgisayar ekranlarda göstermek biraz daha zor olabilir.
library(plotly)
## Zorunlu paket yükleniyor: ggplot2
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
library(plotly)
library(reshape2)
cf.mod <- coef(coklureg)
hsGPA <- seq(min(gpa1$hsGPA),max(gpa1$hsGPA),length.out=50)
ACT <- seq(min(gpa1$ACT),max(gpa1$ACT),length.out=50)
colGPA <- t(outer(hsGPA, ACT, function(x,y) cf.mod[1]+cf.mod[2]*x+cf.mod[3]*y))
p <- plot_ly(x=~hsGPA, y=~ACT, z=~colGPA,type="surface") %>%
add_trace(data=gpa1, x=gpa1$hsGPA, y=gpa1$ACT, z=gpa1$colGPA, mode="markers", type="scatter3d", opacity=0.5)
p
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
stargazer(list(coklureg,basitreg),type = "text")
##
## =================================================================
## Dependent variable:
## ---------------------------------------------
## colGPA
## (1) (2)
## -----------------------------------------------------------------
## hsGPA 0.453***
## (0.096)
##
## ACT 0.009 0.027**
## (0.011) (0.011)
##
## Constant 1.286*** 2.403***
## (0.341) (0.264)
##
## -----------------------------------------------------------------
## Observations 141 141
## R2 0.176 0.043
## Adjusted R2 0.164 0.036
## Residual Std. Error 0.340 (df = 138) 0.366 (df = 139)
## F Statistic 14.781*** (df = 2; 138) 6.207** (df = 1; 139)
## =================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
##Çoklu regresyon modeli şu şekilde yazılabilir
colGPA=β0+β1ACT+β2hsGPA+u
##Çoklu regresyonun sonuçları şu şekilde rapor edilebilir
colGPA=1.286+0.009ACT+0.453hsGPA
sapmareg<-lm(hsGPA~ACT,data = gpa1)
stargazer(list(sapmareg),type = "text")
##
## ===============================================
## Dependent variable:
## ---------------------------
## hsGPA
## -----------------------------------------------
## ACT 0.039***
## (0.009)
##
## Constant 2.463***
## (0.218)
##
## -----------------------------------------------
## Observations 141
## R2 0.120
## Adjusted R2 0.113
## Residual Std. Error 0.301 (df = 139)
## F Statistic 18.879*** (df = 1; 139)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01
standarthata<-summary(coklureg)$sigma
standarthata
## [1] 0.3403158
Standart hatanın karesi bize regresyonumuz hata payının varyansını verecektir.
sigmakare<-standarthata^2
sigmakare
## [1] 0.1158148
Rjkare<-summary(sapmareg)$r.squared
Rjkare
## [1] 0.1195815
Varyans enflasyon faktörü (Variance Inflation Factor (VIF))’i hesaplayalım
VIF<-1/(1-Rjkare)
VIF
## [1] 1.135823
n<-nobs(coklureg)
n
## [1] 141
varxj<-var(gpa1$hsGPA)*(n-1)/n
varxj
## [1] 0.1016267
varbetaj<-(1/n)*(sigmakare/varxj)*VIF
varbetaj
## [1] 0.009180115
sdbetaj<-varbetaj^0.5
sdbetaj
## [1] 0.09581292
VIF’i bulmanın kolay bir yolu car paketini yüklemektir.