As I understand it, I need to first find the residuals and then find estimates for theta but am not sure on how to do either. I am supposed to first assume i.i.d errors and then assume heteroskedasticity but no autocorrelation in the errors.
Is it something of the form:
inv(n)Z’(y-Xtheta))'W(inv(n)Z’(y-Xtheta))
Where W is the weighting matrix for the heteroskedastic case and W=I for i.i.d
Any help would be appreciated,
Thanks
EDIT: Here is what I have so far, not sure how correct it is:
using Optim
using CSV
data=CSV.read(“data Ass_2.csv”, datarow=2)
docvis=Array(data[:docvis])
y=docvis #x and z variables
private=Array(data[:private])
chronic=Array(data[:chronic])
female=Array(data[:female])
income=Array(data[:income])
age=Array(data[:age])
black=Array(data[:black])
hispanic=Array(data[:hispanic])
β0=ones(4412) collection of instruments
Z=[private chronic female age black hispanic]
You define a function f(σ)=y-X̂*θ that does not depend on the input variable σ. So your function is constant. I think you may benefit from a better understanding of how to define methods and functions in Julia, see the resources at Get started with Julia
Optim works well to minimise least squares problems, but you can also consider using LsqFit.jl
By the way, the gmmresults() function is general purpose, for possibly nonlinear in the parameters moment conditions. As such, it’s using Optim.jl, which is not necessary in the case of linear in parameters moment conditions. But, it works.
I’d say GMM is something that doesn’t really need a package. The estimator is defined by two things, the moment conditions, and the weight matrix (or just the moment conditions, if the continuously updating version is used). The following code will compute the ordinary two-step GMM estimator and everything needed to compute the estimated covariance matrix, given a function defining the moment conditions and the weight matrix (fminunc is just a convenience frontend to Optim). There’s not much there to need a dedicated package. I’d say it would be more appropriate to place GMM into a more general econometrics package.
# ordinary GMM with provided weight matrix
# weight should be gXg positive definite
function gmm(moments, theta, weight)
# average moments
m = theta -> vec(mean(moments(theta),dims=1)) # 1Xg
# moment contributions
momentcontrib = theta -> moments(theta) # nXg
# GMM criterion
obj = theta -> ((m(theta))'weight*m(theta))
# do minimization
thetahat, objvalue, converged = fminunc(obj, theta)
# derivative of average moments
D = (Calculus.jacobian(m, vec(thetahat), :central))'
# moment contributions at estimate
ms = momentcontrib(thetahat)
return thetahat, objvalue, D, ms, converged
end
I believe the example by @mcreel is very good. In case you want more of something close that, you may consider my lecture notes/notebook at GMM forth and back
The function gmmresults.jl, mentioned a while ago upthread, uses the gmm code I posted recently, above, to compute the estimator, then it prints out nicely formatted results. To do this, it uses the PrettyTables.jl package. This package allows export to LaTeX. gmmresults.jl doesn’t use this feature, but it should be quite easy to modify the code to add this option.
I’m working on my own implementation of GMM for a project and thought I’d contribute to the many excellent implementations in this thread with a minimal example that makes good use of unicode characters to correspond to familiar mathematical notation. This example assumes you have already programmed a function g(z,θ) that returns an L \times T matrix of L moment contributions computed at T observations (i.e. the t-th column of the output of g is equal to g(z_t,\theta).
using Optim, ForwardDiff
function gmm(g, z, θ₀, W)
# Parameter estimates
f(θ) = (ḡ = mean(g(z,θ), dims=2); dot(ḡ, W, ḡ)) # GMM criterion
result = Optim.optimize(f, θ₀, LBFGS(); autodiff=:forward) # run optimization
θ̂ = Optim.minimizer(result) # access result
# Standard errors: (G'WG)⁻¹G'WΩ̂WG(G'WG)⁻¹ →ᵖ Var(√T̅(θ̂-θ))
T = size(z,1)
G = ForwardDiff.jacobian(θ -> mean(g(z,θ), dims=2), θ̂) # G →ᵖ E[∂g/∂θ′]
D = inv(G'*W*G)
Ω̂ = (ĝ = g(z,θ̂); ĝ*ĝ'/T) # Ω̂ = T⁻¹Σₜ[g(zₜ,θ̂)g(zₜ,θ̂)′] →ᵖ Ω
V̂ = D*G'*W*Ω̂*W*G*D # (G'WG)⁻¹G'WΩ̂WG(G'WG)⁻¹
σ̂ = sqrt.(diag(V̂/T))
return θ̂, σ̂
end
You can test it out by comparing it to OLS using this moment function:
function g(z,θ)
y = @view z[:,1]
X = @view z[:,2:end]
ε̂ = y - X*θ
transpose(X.*ε̂)
end
Understanding the mechanics of actually getting GMM to work was a bit of a hurdle for me coming from Matlab, so I hope this is helpful to others making the switch.
The below is a code for my case
using Pkg
import Pkg;
using DataFrames;
using CSV;
using Lathe;
using GLM;
using Statistics;
using StatsPlots;
using MLBase;
df=DataFrame(CSV.File(“C:Users/GABBYMOSH/OneDrive/Desktop/Work/hsng.csv”)
first(df,5);
names(df) #Multiple linear Model
fm=@formula(rent~hsngval+pcturban)
linearRegressor=lm(fm,df)
##GMM
hsngval=Array(df[:hsngval]);
pcturban=Array(df[:pcturban]);
faminc=Array(df[:faminc]);
reg2=Array(df[:reg2]);
reg4=Array(df[:reg4]); #Instrument variables
Z=[pcturban faminc reg2 reg4];
B0=ones(50);
X=[B0 hsngval pcturban];
rent=Array(df[:rent]);
y=rent;
