Örnek3.1 Kolej Gpa’sının Belirleyicileri ödeviniz olarak da verdiğim kolej Gpa’sının belirleyicileri örneğinin datasını indirelim ve içinde bulunan değişkenlere göz atalım.

library(wooldridge)
data(gpa1)
library(rmarkdown)
paged_table(gpa1)

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. Kitabınızda regresyon sonuçları direk verilen örneğin sonuçlarını bulmak isterseniz şu komutları kullanabilirsiniz.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. Kitabınızda regresyon sonuçları direk verilen örneğin sonuçlarını bulmak isterseniz şu komutları kullanabilirsiniz.

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

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. En azından örneğimizdeki gibi üç boyutlu bir grafiğe ihtiyacımız var. Önce bütün değişkenler için üç boyutlu bir dağılım grafiği çizelim daha sonra doğru yerine bir düzlem eklememiz gerekecek.

library(plotly)
## Loading required package: 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(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

Eğitim yetenek arttığı için artacağından, yeteneği modele katmamak

β1~ tahmininin o kadar artmasına yol açar. Şimdi bu hesaplamamızı bir örnek üzerinden gösterelim. Bu örneği hatırlarsak

##install.packages("stargazer")
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

δ1~iki bağımsız değişkenin birbirleri arasındaki basit regresyon ilişkisinden de bulunabilir

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

Örneğimizi kullanarak bu parçaları çıkaralım.

standarthata<-summary(coklureg)$sigma
standarthata
## [1] 0.3403158

Standart hatayi hesaplamak için sadece sigma komutunu kullanıyoruz. Standart hatayi colGPA örneği için 0.34 olarak hesapladık. σ2 yani standart hatanın karesi bize regresyonumuzun hata payının varyansını verecektir

sigmakare<-standarthata^2
sigmakare
## [1] 0.1158148

Regresörlerimiz arasındaki Rj2 hesaplamak için daha önce hesapladığımız sapmareg’i kullanabiliriz. Bu Rj2bizim için lise not ortalamasıyla ACT sonuçları arasındaki ilişkiyi gösterecektir. Eğer Rj2 bir’e yakınsa çoklu eş doğrusallık muhtemel olabilir.

Rjkare<-summary(sapmareg)$r.squared
Rjkare
## [1] 0.1195815

İki bağımsız değişkenimizin basit regresyonundan elde ettiğimiz Rj2yaklaşık olarak yüzde 12’dir.11−Rj2varyans enflasyon faktörü (Variance inflation factor (VIF))’i hesaplayalım.

VIF<-1/(1-Rjkare)
VIF
## [1] 1.135823

son olarak n’i bulmalıyız.

n<-nobs(coklureg)
n
## [1] 141

son olarak var(xj)’yi bulmamız lazım.

varxj<-var(gpa1$hsGPA)*(n-1)/n
varxj
## [1] 0.1016267

Şimdi formülümüzü kullanarak var(βj^)’yi hesaplayabiliriz

varbetaj<-(1/n)*(sigmakare/varxj)*VIF
varbetaj
## [1] 0.009180115

βj’nin standart hatası var(βj^)’nin karaköküdür.

sdbetaj<-varbetaj^0.5
sdbetaj
## [1] 0.09581292

Tablodan hsGPA’in standart hatasına bakarsanız bulduğumuz rakamla aynı olduğunu göreceksiniz. VIF’i bulmanın kolay bir yolu car paketini yüklemektir.

require(car)
## Loading required package: car
## Loading required package: carData
vif(coklureg)
##    hsGPA      ACT 
## 1.135823 1.135823