FGLS and Other Estimators

The feasible GLS estimator of the AR(p) model can be estimated using gretl in a number of ways. For first order autocorrelated models the ar1 command can be used. There are a num­
ber of estimators available by option including the Cochrane-Orcutt (iterated), the Prais-Winsten (iterated), and the Hildreth-Lu search procedure. Examples are:

1 list x = d_u const

2 ar1 inf x # Cochrane-Orcutt (default)

3 ar1 inf x —pwe # Prais-Winsten

4 ar1 inf x —hilu —no-corc # Hildreth-Lu

The results are collected in a model table below.

AR(1) Errors
Dependent variable: inf

(CO)

(PW)

(HL)

const

0.7609**

0.7862**

0.7608**

(0.1238)

(0.1218)

(0.1245)

d_u

-0.6944**

-0.7024**

-0.6953*

(0.2429)

(0.2430)

(0.2430)

P

0.55739

0.55825

.56

n

89

90

89

R2

0.3407

0.3418

0.3406

Standard errors in parentheses
* indicates significance at the 10 percent level
** indicates significance at the 5 percent level
CO = Cochrane Orcutt, PW=Prais-Winsten, HL=Hildreth-Lu

You can see that there are minor differences produced by these options. If the —no-corc option is not used with –hilu then the Hildreth-Lu estimator is modified slightly to perform additional iterations as the end. Notice that the Prais-Winsten is the only procedure to use all 90 observations.

For higher order models there are two commands worth taking note of. The ar command estimates a linear regression with arbitrary autocorrelation structure. It uses a generalization of the Cochrane-Orcutt iterative procedure to obtain estimates.

The other estimator is arima, the syntax for which appears below:

Подпись: Arguments: Options: Подпись: Examples:

Подпись: arima

p d q [ ; P D Q ] ; depvar [ indepvars ]

—verbose [print details of iterations)

—vcv (print covariance matrix)

–nessian (see below)

—opg (see below)

—nc (do not include a constant)

—conditional (use conditional maximum likelihood) —x-12-arima (use X-12-ARIMA for estimation)

— lbf gs (use L-BFGS-B maximizer)

—y-dif f-only (ARIMAX special, see below) –save-enat (see below) arima 1 0 2 ; у

arima 202 ; у 0 xl x2 –verbose arima Oil; Oil; у —nc

The default estimation method for arima in gretl is to estimate the parameters of the model using the “native” gretl ARMA functionality, with estimation by exact maximum likelihood using the Kalman filter. You can estimate the parameters via conditional maximum likelihood as well.

Estimating the simple AR(1) regression using these estimators is done:

1 ar 1 ; inf x

2 arima 100; inf x

For the ar command, list the lag numbers for the desired residuals. In the case of AR(1) this is just 1. This is followed by a semicolon and then the regression to estimate. The arima syntax is similar, except you specify p, d, and q, where p is the order of the desired autocorrelation, d is the number of differences to take of the time-series, and q is the order of any moving average terms you might have in the residuals.

The outcome for the simple ARIMA(1,0,0) ia

ARMAX, using observations 1987:2-2009:3 (T = 90)

Estimated using Kalman filter (exact ML)

Dependent variable: inf Standard errors based on Hessian

coefficient std. error z p-value

const

0.786212

0.120601

6.519

7.07e-011

***

phi_1

0.558827

0.0877359

6.369

1.90e-010

***

d_u

-0.702558

0.242234

-2.900

0.0037

***

S. D. dependent var

0.636819

S. D. of innovations

0.510937

Akaike criterion

142.9118

Hannan-Quinn

146.9441

maginary Modulus

Frequency

0.0000 1.7895

0.0000

Mean dependent var 0.791111 Mean of innovations -0.003996 Log-likelihood -67.45590

Schwarz criterion 152.9110

Real

AR

Root 1 1.7895

These are very similar to the ones above. The coefficient labeled phi_1 is the estimate of the autocorrelation parameter. The root of this equation is 1/phi_1. The roots (or modulus) must be greater than 1 in absolute value in order for the model to be stationary.

9.11 Script

1 open "@gretldirdatapoeokun. gdt"

2 set echo off

3 # change variable attributes

4 setinfo g – d "percentage change in U. S. Gross Domestic Product, seasonally

5 adjusted" -n "Real GDP growth"

6 setinfo u – d "U. S. Civilian Unemployment Rate (Seasonally adjusted)" – n

7 "Unemployment Rate"

8

8 # plot series and save output to files

9 gnuplot g —with-lines —time-series —output="@workdirokun_g. plt"

10 gnuplot u —with-lines —time-series —output="@workdirokun_u. plt"

12

11 # graphing multiple time-series

12 scatters g u –with-lines

15

13 diff u

14 setinfo d_u – d "Change in U. S. Civilian Unemployment

15 Rate (Seasonally adjusted)" – n

16 "D. Unemployment Rate"

17 scatters g d_u –with-lines –output=display

21

18 # distributed lag models

19 ols d_u const g(0 to -3)

20 smpl 1986:1 2009:3

21 ols d_u const g(0 to -2)

26

22 gnuplot g g_1

28

23 # correlogram and confidence interval

24 corrgm g 12

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

matrix ac = corrgm(g, 12) matrix lb = ac[,1]-1.96/sqrt($nobs) matrix ub = ac[,1]+1.96/sqrt($nobs) matrix all = lb~ac[,1]~ub colnames(all, "Lower AC Upper ")

printf "nAutocorrelations and 95%% confidence intervalsn %9.4fn", all

# Phillips curve

open "@gretldirdatapoephillips_aus. gdt" diff u

setinfo inf – d "Australian Inflation Rate" – n "Inflation Rate" setinfo d_u – d "Change in Australian Civilian

Unemployment Rate (Seasonally adjusted)" – n

"D. Unemployment Rate" scatters inf d_u –with-lines

ols inf const d_u series ehat = $uhat gnuplot ehat —time-series corrgm ehat

# LM tests

ols ehat const d_u ehat(-1) scalar NR2 = $trsq pvalue X 1 NR2

ols ehat const d_u ehat(-1 to -4) scalar NR2 = $trsq pvalue X 4 NR2

ols inf const d_u modtest 1 —autocorr modtest 4 —autocorr —quiet

# HAC standard errors

open "@gretldirdatapoephillips_aus. gdt" set hac_kernel bartlett set hac_lag nw2 diff u

ols inf const d_u modeltab add

ols inf const d_u –robust modeltab add modeltab show modeltab free

# nonlinear least squares estimation of regression w/AR(1) errors open "@gretldirdatapoephillips_aus. gdt"

diff u

ols inf const d_u –quiet

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

scalar betal = $coeff(const) scalar beta2 = $coeff(d_u) scalar rho = 0

nls inf = beta1*(1-rho) + rho*inf(-1) + beta2*(d_u-rho*d_u(-1)) params rho beta1 beta2 end nls

scalar delta = $coeff(beta1)*(1-$coeff(rho)) scalar delta1 = -$coeff(rho)*$coeff(beta2)

printf "nThe estimated delta is %.3f and the estimated delta1 is %.3f.n",delta, delta1 scalar sser=$ess

# estimation of more general model ols inf const inf(-1) d_u(0 to -1) scalar sseu=$ess

scalar fstat = (sser-sseu)/(sseu/$df) pvalue X 1 fstat pvalue F 1 $df fstat omit d_u(-1)

ols inf const inf(-1) d_u(0 to -1) modeltab add

ols inf const inf(-1) d_u(0) modeltab add modeltab show modeltab free

# model selection function

function matrix modelsel (series y, list xvars) ols y xvars —quiet scalar sse = $ess scalar N = $nobs scalar K = nelem(xvars) scalar aic = ln(sse/N)+2*K/N scalar bic = ln(sse/N)+K*ln(N)/N matrix A = { K, N, aic, bic} printf "nRegressors: %sn",varname(xvars)

printf "K = %d, N = %d, AIC = %.4f SC = %.4f.n",K, N,aic, bic return A end function

# using the modelsel function

list x = const inf(-1) d_u(0 to -1)

matrix a = modelsel(inf, x)

list x0 = const

matrix b = modelsel(inf, x)

list x = const d_u inf(-1)

# putting the model selection results into a matrix open "@gretldirdatapoephillips_aus. gdt"

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

diff u

smpl 1988:3 2009:3 matrix A = {} scalar q = 0 loop p = 1..6 —quiet if p = 1

list x = const inf(-1) d_u

else

list x = const inf(-1 to -$p) d_u endif

matrix a = $p~q~modelsel(inf, x) matrix A = A | a modelsel(inf, x) endloop scalar q = 1 loop p = 1..6 —quiet if p = 1

list x = const inf(-1) d_u(0 to -1)

else

list x = const inf(-1 to -$p) d_u(0 to -1) endif

matrix a = $p~q~modelsel(inf, x) matrix A = A | a endloop

colnames(A,"p q K N AIC SC ") print A

smpl full

ols inf const inf(-1 to -4) d_u —robust

# improved modelsel2 function for ARDL function matrix modelsel2 (list xvars)

ols xvars –quiet

scalar sse = $ess

scalar N = $nobs

scalar K = nelem(xvars)-1

scalar aic = ln(sse/N)+2*K/N

scalar bic = ln(sse/N)+K*ln(N)/N

matrix A = { K, N, aic, bic}

# printf "nDependent variable and Regressors: %sn",varname(xvars)

# printf "K = %d, N = %d, AIC = %.4f SC = %.4f.n",K, N,aic, bic return A

end function

# using modelsel2

open "@gretldirdatapoeokun. gdt" diff u

smpl 1986:1 2009:3

matrix A = {}

loop p = 0..2 —quiet

loop q = 1..3 —quiet

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

217

218

219

220

221

222

223

224

225

226

227

228

229

230

231

232

233

234

if p=0

list vars = d_u g(0 to – q) const

else

list vars = d_u(0 to – p) g(0 to – q) const endif

matrix a = p~q~modelsel2(vars) matrix A = A | a endloop endloop

colnames(A,"p q K N AIC SC ") print A

function modelsel clear smpl full

ols d_u(0 to -1) g(0 to -1) const loop i=1..4

modtest $i —autocorr —quiet endloop

open "@gretldirdatapoeokun. gdt" smpl 1986:3 2009:3 matrix A = {} scalar q=0

loop p = 1..5 —quiet

list vars = g(0 to – p) const matrix a = p~q~modelsel2(vars) matrix A = A | a endloop

colnames(A,"p q K N AIC SC ") print A

function modelsel clear

# loop to test for autocorrelation in ARDL open "@gretldirdatapoephillips_aus. gdt" diff u

ols inf(0 to -1) d_u const loop i=1..5

modtest $i —autocorr —quiet endloop

# loop to test for autocorrelation at several lags open "@gretldirdatapoeokun. gdt"

ols g(0 to -2) const series res = $uhat corrgm res loop i = 1..4

modtest $i —autocorr —quiet endloop

# model selection for Okun data open "@gretldirdatapoeokun. gdt"

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255

256

257

258

259

260

261

262

263

264

265

266

267

268

269

270

271

272

273

274

275

276

277

278

279

280

281

282

283

284

285

smpl 1986:3 2009:3 matrix A = {} scalar q=0

loop p = 1..5 —quiet

list vars = g(0 to – p) const matrix a = p~q~modelsel2(vars) matrix A = A | a endloop

colnames(A,"p q K N AIC SC ") print A

# estimation of preferred model and a forecast open "@gretldirdatapoeokun. gdt"

ols g(0 to -2) const dataset addobs 3

fcast 2009:4 2010:2 —plot="@workdirar2plot1.plt"

# multiplier analysis

open "@gretldirdatapoeokun. gdt"

matrix y = { g }

scalar T = $nobs

matrix sm1 = zeros(T,1)

scalar a = .38

smpl 1 round((T+1)/2)

scalar stv = mean(y)

smpl full

loop i=1..T –quiet if i = 1

matrix sm1[i]=stv

else

matrix sm1[$i]=a*y[$i] + (1-a)*sm1[i-1] endif endloop

series exsm = sm1

gnuplot g exsm —time-series

scalar a = .8 loop i=1..T –quiet if i = 1

matrix sm1[i]=stv

else

matrix sm1[$i]=a*y[$i] + (1-a)*sm1[i-1] endif endloop

series exsm8 = sm1

gnuplot g exsm8 —time-series

open "@gretldirdatapoeokun. gdt" diff u

ols d_u(0 to -1) g(0 to -1) const scalar b0 = $coeff(g)

286

287

288

289

290

291

292

293

294

295

296

297

298

299

300

301

302

303

304

305

306

307

308

309

310

311

312

313

314

315

316

317

318

319

320

321

322

323

324

325

326

327

328

329

330

331

332

333

334

335

336

scalar b1 =$coeff(d_u_1)*b0+$coeff(g_1) scalar b2 = b1*$coeff(d_u_1) scalar b3 = b2*$coeff(d_u_1)

# Matrix & Series Plot

open "@gretldirdatapoeokun. gdt" diff u

ols d_u(0 to -1) g(0 to -1) const scalar h = 8

matrix mult = zeros(h,2) loop i=1..h

mult[i,1] = i-1 scalar b0 = $coeff(g)

scalar b1 = $coeff(d_u_1)*b0+$coeff(g_1) if i=1

mult[i,2]=b0 elif i=2

mult[i,2]=b1

else

mult[i,2]=mult[i-1,2]*$coeff(d_u_1)

endif

endloop

gnuplot 2 1 —matrix=mult —output=display —with-lines —suppress-fitted

printf "nThe impact and delay multipliers are n %10.5fn", mult

nulldata 8 –preserve series m = mult[,2] series lag = mult[,1]

setinfo m – d "Multipliers" – n "Multiplier"

gnuplot m index —with-lines —output=display —suppress-fitted

# appendix

open "@gretldirdatapoephillips_aus. gdt" diff u

setinfo inf – d "Australian Inflation Rate" – n "Inflation Rate" setinfo d_u – d "Change in Australian Civilian

Unemployment Rate (Seasonally adjusted)" – n

"D. Unemployment Rate"

# Durbin-Watson with p-value list x = d_u const

ols inf x

scalar dw_p = $dwpval print dw_p

# various ways to estimate AR(1) regression ar1 inf x

modeltab add ar1 inf x –pwe

337 modeltab add

338 ar1 inf x —hilu —no-corc

339 modeltab add

340 modeltab show

341 modeltab free

342

342 ar 1 ; inf x

343 arima 10 0; inf x

residual

0 2 4 6 8 to 12

lag

Figure 9.11: This plot shows that the residuals from the simple Phillips curve model are serially correlated. Australia, 1987:1 – 2009:3.

Breusch-Godfrey teat for autocorrelation up to order 4 OLS, using observations 1987:2-2009:3 (T = 90) Dependent variable: uhat

COn3t

d_u uhat_l uhat_2 uhat_3 l. lll 4

Unadjusted R-squared = 0.407466

Test statistic: LMF – 14.440976,

with P-value = P(F(4,84) > 14.441) = S. lSe-009

Alternative statistic: TR~2 = 36.671897,

with p-value = P(Chi-square(4) > 36.6719) = 2.ІЄ-007

Ljung-Box Q1 = 82.4327,

with p-value = P(Chi-square(4) > 82.4327) = 5.3ІЄ-017

Figure 9.12: Using Test>Autocorrelation from the model pull-down menu will generate the following output. The alternative hypothesis is AR(4).

Using numerical derivatives

Tolerance = 1.81899e-012

Convergence achieved after 21 iterations

Model 2: NL5, using observations 1987:3-2009:3 (T = 89) inf = betal*(1-rho) + rho*inf(-l) + beta2*<d_u-rho*d_u(-1))

estimate std. error t-ratio p-value

Figure 9.13: Nonlinear least squares results for the AR(1) regression model.

Figure 9.14: Residual correlogram for Okun AR(2)

Figure 9.15: Using Data>Add observations from the main gretl pull-down menu will extend the sample period. This is necessary to generate forecasts.

Figure 9.16: Forecast dialog box

Figure 9.19: Impact and delay multipliers for an ARDL(1,1) of the change in unemployment caused by 1% increase in U. S. GDP growth.

Leave a reply

You may use these HTML tags and attributes: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>