Baltagi Table 8.1 Replication with plm I saw this previous post from 2017 that talks about the same table but only focuses on OLS.  I'm also working on replicating Table 8.1 in Baltagi's Economic Analysis of Panel data book.  I can replicate the table in Stata and I can replicate most of the table using the plm package in R. Unfortunately, I cannot replicate the Arellano-Bond Twostep GMM values.  I can confirm that I can get the same results as the Table using Stata but do not get the same results as those in the table using plm.  I am confused what option or technique is causing the issue.
I've included the R code here for both the one-step and two-step GMM estimates.  As shown, when estimating the one-step, the results match while when estimating the two-step they do not.  I'm also including Table 8.1 at the bottom to show the results presented in the book.
library(plm)
#Loading the Baltagi dataset provided as part of the plm package
data(Cigar)

#transforming to real variables as required
Cigar$real_c <- with(Cigar, (sales*pop)/pop16)
Cigar$real_p <- with(Cigar, (price/cpi)*100)
Cigar$real_pimin <- with(Cigar, (pimin/cpi)*100)
Cigar$real_ndi <- with(Cigar, (ndi/cpi))



gmm_onestep <- pgmm(log(real_c) ~ lag(log(real_c), 1) + log(real_p) + log(real_pimin) + log(real_ndi) | 
                      lag(log(real_c), 2:99), data = Cigar, effect = "twoways", model = "onestep")

#match w/ Table 8.1: 0.84, -0.377, -0.016, 0.14
round(summary(gmm_onestep)$coef, 3)

#Two step does not match Table 8.1: 0.80, -0.379; -0.020, 0.24
gmm_twostep <- pgmm(log(real_c) ~ lag(log(real_c), 1) + log(real_p) + log(real_pimin) + log(real_ndi)  | 
                      lag(log(real_c), 2:99), data = Cigar, effect = "twoways", model = "twostep")

#pgmm reports 0.632, -0.358, -0.002, 0.386
round(summary(gmm_twostep)$coef, 3)


Any help would be greatly appreciated!

 A: The key to obtain similar results is to notice this part in stata code dum3-dum29 in command:
# stata
xtset state year
xtabond lnc lnrp lnrpn lnrdi dum3-dum29, lag(1) twostep 

Hence, we have to model time dummies manually (i.e. use the packege fastDummies).
Baltagi in Stata drops not only two first periods of time dummies, but also the last one (I've commented out them in the code below). You may uncomment two lines to get you initial results with automatic twoways effects. As far as we model time manually, I've changed effect option to individual.
library(fastDummies)
Cigar <- dummy_cols(.data = Cigar, select_columns = "year")

gmm_tmp <- pgmm(log(real_c) ~ lag(log(real_c), 1) + log(real_p) + log(real_pimin) + log(real_ndi) 
                # + year_63 + year_64 
                + year_65 + year_66 + year_67 + year_68 + year_69 + year_70 + year_71 + year_72 +
                  year_73 + year_74 + year_75 + year_76 + year_77 + year_78 + year_79 +
                  year_80 + year_81 + year_82 + year_83 + year_84 + year_85 + year_86 +
                  year_87 + year_88 + year_89 + year_90 + year_91 
                # + year_92
                | lag(log(real_c), 2:99)
                , data = Cigar, effect = "individual", model = "twosteps", transformation = "d"
                )

round(coefficients(gmm_tmp)$coef[1:4], 2)
# 0.80               -0.37               -0.08                0.18

Below is a part of summary output to see that only first two coefficient are significant, indeed (so the difference of two last may stem from numeric differences in matrix procedures).
# Coefficients:
#                          Estimate Std.Error z-value Pr(>|z|)    
#   lag(log(C_real), 1)  0.7991766  0.1962817  4.0716 4.670e-05 ***
#   log(P_real)         -0.3696419  0.0921545 -4.0111 6.043e-05 ***
#   log(Pn_real)        -0.0782846  0.1524425 -0.5135  0.607577    
#   log(Y_real)          0.1846966  0.1215665  1.5193  0.128686  

Notice also that transformation = "d" corresponds to xtabond in Stata, while xtabond2 counterpart seems to be transformation = "ld" in R.
Hope that helps.
