22NisanDers
R Markdown
library(readxl)
tarimdataset <- read_excel("E:/YuksekLisans/Alternatif Regresyon/tarimdataset1.xlsx")
library(DT)
datatable(tarimdataset)Veri setini incelediğimizde 4 değişkene sahip olduğu görülmektedir.
KTKDGR: (cari $), Tarımsal Katma Değer Karekökü
MKN: (adet) , Tarımda kullanılan tarımsal makinelerin, paletli ve tekerlekli traktörlerin sayısıdır.
İTH: (Mal İthalatının %’si), Tarımsal hammadde ithalatının toplam mal ithalatı içersindeki yüzdesidir
GÜB: (ton) Kilometre kare başına gübre tüketimi
Tanımlayıcı istatistiklere göz atalım
summary(tarimdataset)## ktkdgr mkn ithh gubb
## Min. : 4200 Min. : 0.4 Min. : 492.5 Min. : 0.2
## 1st Qu.: 28950 1st Qu.: 37.0 1st Qu.:1088.8 1st Qu.: 227.3
## Median : 47700 Median : 300.9 Median :1558.6 Median : 1566.3
## Mean : 65120 Mean : 2073.0 Mean :1951.5 Mean : 6819.4
## 3rd Qu.: 90950 3rd Qu.: 1656.1 3rd Qu.:2317.1 3rd Qu.: 9553.6
## Max. :240000 Max. :20280.0 Max. :7420.2 Max. :67730.0Çoklu doğrusal regresyon uygulanmadan önce veri setinin yeterli olup olmadığını sorgulamak gerekmektedir. Açıklayıcı değişkenler arası korelasyon tablosuna bir göz atalım:
cor(tarimdataset)## ktkdgr mkn ithh gubb
## ktkdgr 1.0000000 0.7425557 0.4249229 0.7348492
## mkn 0.7425557 1.0000000 0.1811087 0.5459359
## ithh 0.4249229 0.1811087 1.0000000 0.3013769
## gubb 0.7348492 0.5459359 0.3013769 1.0000000library(corrplot)## corrplot 0.84 loadedcorrplot(cor(tarimdataset), type = "upper", order = "hclust",
tl.col = "black", tl.srt = 45)
Açıklayıcı değişkenler arasında yüksek korelasyonlar mevcut fakat bunları şu aşamada göz ardı edeceğiz. Bu arada 0.8’den daha yüksek korelasyon olmadığı da gözlerden kaçmamaktadır.
model<- lm(ktkdgr ~ mkn + ithh + gubb, data = tarimdataset)
summary(model)##
## Call:
## lm(formula = ktkdgr ~ mkn + ithh + gubb, data = tarimdataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -66921 -14990 -3399 14106 64891
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.488e+04 6.107e+03 4.074 0.00015 ***
## mkn 5.797e+00 9.733e-01 5.956 1.90e-07 ***
## ithh 8.149e+00 2.691e+00 3.028 0.00374 **
## gubb 1.807e+00 3.698e-01 4.887 9.24e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 26120 on 55 degrees of freedom
## Multiple R-squared: 0.748, Adjusted R-squared: 0.7343
## F-statistic: 54.43 on 3 and 55 DF, p-value: < 2.2e-16Modelin F istatistiği için hesaplanan p değerine göre kurulan model bütünüyle anlamlıdır. Katsayıların anlamlılığını da tek tek inceleyecek olursak t değerleri -+ 1,96 arasında bulunmadığı için sıfırdan farklıdır denir. Yani katsayılar bireysel olarak da modele bir katkı sağlamaktadır. Bunu bir de regresyon katsayılarının güven aralığı tablosuyla görelim:
guven_araliklari <- confint(model, level = 0.95)
datatable(guven_araliklari)Residual Standard error değerini elle hesaplayalım:
prediction <- predict(model, interval = "prediction")## Warning in predict.lm(model, interval = "prediction"): predictions on current data refer to _future_ responsesprediction <- datatable(prediction)
err <- sqrt(sum((prediction$fit- tarimdataset$ktkdgr)^2)/(59-4))
show(err)## [1] 0(r’da bu değer 26120.5223761 olarak çıkıyor. Buradaki 0 değerine bakmamak lazım) Hataların ortalamasının 0 olup olmadığını inceleyelim.
mean_error_terms <- mean(model$residuals)
mean_error_terms## [1] -1.188774e-12Hataların ortalaması 0’a çok yakındır. Bu varsayım da sağlanmıştır denilebilir. Hataların normal dağılıp dağılmadığını inceleyelim.
library(data.table)
error_terms <- (prediction$fit - tarimdataset$ktkdgr)
error_terms <- data.table(error_terms)
shapiro.test(model$residuals)##
## Shapiro-Wilk normality test
##
## data: model$residuals
## W = 0.96859, p-value = 0.1307Ho hipotezi reddedilemez yani normallik sağlanmış olur.
Multicollinearity’nin problem olup olmadığını VIF ile test edeceğiz. Korelasyon matrisine göre önemli bir problem çıkmamıştı, fakat VIF ile de test edelim. (Bu kod R’da çalışıyor, html formatında burada yazılmıyor) vif(model)
Anova tablosunu oluşturalım
anova(model)## Analysis of Variance Table
##
## Response: ktkdgr
## Df Sum Sq Mean Sq F value Pr(>F)
## mkn 1 8.2118e+10 8.2118e+10 120.359 1.820e-15 ***
## ithh 1 1.2989e+10 1.2989e+10 19.038 5.687e-05 ***
## gubb 1 1.6297e+10 1.6297e+10 23.886 9.238e-06 ***
## Residuals 55 3.7525e+10 6.8228e+08
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Bağımlı değişken ile bağımsız değişkenler arasında lineerlik konusunda bir problem olup olmadığını inceleyeceğiz.
plot(tarimdataset)
Yeni değerler için tahmin yapma ve yüzde 95 güven aralığını belirle
newdata<-data.frame(mkn=2000, ithh=2000, gubb=7000)
aralik_tahmini <- predict(model, newdata, interval = "confidence")
aralik_tahmini## fit lwr upr
## 1 65418.16 58595.07 72241.25library(readxl)
MaasDataset <- read_excel("E:/YuksekLisans/LisansDersi/MaasDataset.xlsx")
library(DT)
datatable(MaasDataset)Dummy değişkenini manuel olarak oluşturalım. R bunu otomatik olarak yapıyor fakat bu daha sonra gösterilecektir. BOLGE değişkeninin 1 değerini aldığı Dummy D1 olarak, 2 değerini aldığı Dummy değişkeni D2 olarak kodlanmıştır. Yalın hali ise BOLGE değişkeninin 3 değerini aldığı durumlar olmuş oluyor. Bunu R’a tanıtıp regresyon modelini oluşturup özetini R’dan isteyelim.
D1 <- (MaasDataset$BOLGE == "1")*1
D2 <- (MaasDataset$BOLGE=="2")*1
model_maas <- lm(ODEME ~ HARCAMA + D1+D2, data=MaasDataset)
summary(model_maas)##
## Call:
## lm(formula = ODEME ~ HARCAMA + D1 + D2, data = MaasDataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3937.9 -1517.5 -158.2 1291.6 6133.2
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13265.2424 1394.9433 9.510 1.58e-12 ***
## HARCAMA 3.2898 0.3176 10.357 1.02e-13 ***
## D1 -1672.5561 800.9841 -2.088 0.0422 *
## D2 -1143.4021 860.9266 -1.328 0.1906
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2270 on 47 degrees of freedom
## Multiple R-squared: 0.7228, Adjusted R-squared: 0.7051
## F-statistic: 40.85 on 3 and 47 DF, p-value: 3.833e-13library(car)## Loading required package: carDataAnova(model_maas)## Anova Table (Type II tests)
##
## Response: ODEME
## Sum Sq Df F value Pr(>F)
## HARCAMA 552604863 1 107.2772 1.018e-13 ***
## D1 22460579 1 4.3603 0.04223 *
## D2 9086006 1 1.7639 0.19056
## Residuals 242105779 47
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Bölgelere göre X’e karşı Y’nin grafiğini çizdirelim
dummy <- factor(MaasDataset$BOLGE)
plot(MaasDataset$HARCAMA, MaasDataset$ODEME, type = "n")
text(MaasDataset$HARCAMA, MaasDataset$ODEME ,substring(as.character(dummy),1,1))
Uyum değerlerine karşın artıklar ve normal olasılık grafiğini çizdirelim.
layout(matrix(c(1,2),1,1))
plot(fitted(model_maas), resid(model_maas),type = "n")
text(fitted(model_maas), resid(model_maas), substring(as.character(MaasDataset$BOLGE),1,1))
abline(0,0)
qqnorm(resid(model_maas))
qqline(resid(model_maas))