library(WDI)
df = WDI(indicator=c(aay='AG.LND.ARBL.ZS', kbh='AG.LND.ARBL.HA.PC' ), country=c('TR'),  start=1999, end=2020)
library(dynlm)
## Warning: package 'dynlm' was built under R version 4.1.3
## Zorunlu paket yükleniyor: zoo
## Warning: package 'zoo' was built under R version 4.1.3
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
 df.ts <- ts(df, start=c(1999), end=c(2020),frequency=1)
 df.ts<-df.ts[,c("aay","kbh")]
 plot(df.ts[,"aay"], ylab="ekilebilir.arazi.alanının.yüzdesi")

 plot(df.ts[,"kbh"], ylab="ekilebilir.arazi.kişi.başı.hektar")

gecikmeliaay <- data.frame(cbind(df.ts[,"aay"], lag(df.ts[,"aay"],-1 )))
names(gecikmeliaay) <- c("aay","aay1gec")
head(gecikmeliaay)
##        aay  aay1gec
## 1 31.54633       NA
## 2 30.95773 31.54633
## 3 30.92525 30.95773
## 4 31.17602 30.92525
## 5 30.36784 31.17602
## 6 31.01620 30.36784
ekilebilirarazialanınınyüzdesi <- df.ts[,"aay"]
deltaekilebilirarazialanınınyüzdesi <- diff(df.ts[,"aay"])
plot(ekilebilirarazialanınınyüzdesi)

plot(deltaekilebilirarazialanınınyüzdesi)

phillips.reg <- dynlm(ekilebilirarazialanınınyüzdesi~deltaekilebilirarazialanınınyüzdesi)
summary(phillips.reg)
## 
## Time series regression with "ts" data:
## Start = 2000, End = 2018
## 
## Call:
## dynlm(formula = ekilebilirarazialanınınyüzdesi ~ deltaekilebilirarazialanınınyüzdesi)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.7390 -1.7456 -0.6646  2.2035  2.6946 
## 
## Coefficients:
##                                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                          28.5240     0.5520   51.68   <2e-16 ***
## deltaekilebilirarazialanınınyüzdesi   0.4432     0.9429    0.47    0.644    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.037 on 17 degrees of freedom
##   (0 observations deleted due to missingness)
## Multiple R-squared:  0.01283,    Adjusted R-squared:  -0.04524 
## F-statistic: 0.2209 on 1 and 17 DF,  p-value: 0.6443
uhat <- resid(phillips.reg)
plot(uhat)

acf(uhat)

library(lmtest)
## Warning: package 'lmtest' was built under R version 4.1.3
bgtest(phillips.reg, order=1, type="F")
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  phillips.reg
## LM test = 94.307, df1 = 1, df2 = 16, p-value = 4.123e-08