Fit measures for GMM Arellano-Bond estimator in R A colleague and I have been working with difference GMM, i.e. the Arellano-Bond estimator, in R. Our option has been to use the pgmm command from the plm package. However, now I am struggling to test the fit of my models, since the package plm itself does not bring such functionality.
Would there be a good way to test for the model fit of a difference-GMM model? Of course, an approach would be to test the square of the correlation between Y and Yhat, as in Bloom, Bond et al. (2001). However, the pgmm command does not give the predictions, only the fitted-values.
Could anyone give any tips on that?
 A: If you are using the pgmm command from the package plm you do have good ways to assess model fit available to you (if this is not what you're using it's what you should be using). The standard way to test model fit would be to look at the $J$-test of overidentifying restrictions, also called the Sargan test (see here or here). Strictly speaking this is a test of the validity of your all instruments (given that some subset of your instruments is valid). This test is listed under "Sargan test" when you use the summary command on your pgmm model. summary also shows test of serial correlation of the residuals.
I gather what you want is less to test the Arellano-Bond model itself and more to test the degree to which the predicted time series seems to be close to the result. You can actually reproduce the square correlation used in Bloom et al. (2001) pretty easily. The fitted values are exactly what you need. Starting with an example from the plm documentation found here:
data("EmplUK", package = "plm")
## Arellano and Bond (1991), table 4b 
z1 <- pgmm(log(emp) ~ lag(log(emp), 1:2) + lag(log(wage), 0:1)
           + log(capital) + lag(log(output), 0:1) | lag(log(emp), 2:99),
            data = EmplUK, effect = "twoways", model = "twosteps")

# Shows Sargan test and serial correlation of residuals.
summary(z1)

# Getting the actual Y values out. There is likely a more elegant way to do this.
Y <- c()
for( i in 1:length(z1$model)){ Y <- c(Y,z1$model[[i]][,1])}

# Note that the fitted-values are exactly what you need.
Yhat <- z1$fitted.values[1:length(z1$fitted.values)]

# Squared correlation of fitted values and actual values
cor(Y,Yhat)^2

