Estimación de la productividad total de los factores (PTF) en R

Question

Estoy tratando de calcular la PTF con el estprod en R porque este paquete me permite calcular con Gross-Output (como con levinsohn petrin) pero tengo problemas para conseguir omega

Mientras que el prodest paquete puede calcular omega pero no puedo sumar Producción bruta .

¿Alguien ha calculado la PTF en R?

[actualizado] Ejemplo reproducible

library(prodest)
library(estprod)

data = structure(list(Y = c(12.8268952320028, 12.8521337216858, 12.8297798993317, 12.6357098405348, 
                            12.4004643557221, 12.3600068590732, 12.826153741752, 12.8990697648315, 
                            12.8393196511736, 12.7619091581355), 
                      Labor = c(18.2268433780233, 18.1280743937883, 
                            16.2004953940248, 18.4010948157752, 16.1374389570964, 18.3902253199657, 
                            16.2359220634517, 16.2468898407345, 16.2720069307031, 16.2673373011924), 
                      Capital = c(12.1642906069818, 12.088450259469, 
                           12.0237634086603, 12.0750443461238, 12.0098651639831, 12.0601371534517, 
                           11.9679864443247, 11.9299321784366, 12.0016971178669, 12.1209493577051), 
                      Materials = c(18.1472736836291, 17.7735308597061, 
                          17.1908332586684, 18.3914858811768, 18.2747180360049, 16.9304018090212, 
                          18.4414110290762, 18.0188590654862, 18.1427756804037, 15.805209131618), 
                      id = c(1070101, 1070101, 1070101, 1070101, 
                              1070101, 1070101, 1070102, 1070102, 1070102, 1070102), 
                      year = c(2013,2014, 2015, 2016, 2017, 2018, 2013, 2014, 2015, 2016)), 
                 class = c("tbl_df", "tbl", "data.frame"), row.names = c(NA, -10L))

# Y, Labor, Capital, Materials are in logarithms
mod1 = estprod::levinsohn_petrin(data,
                                 formula = Y ~ Labor | Capital | Materials, 
                                 id = "id", 
                                 time = "year", 
                                 reps = 20,
                                 gross = T) #estprod with Gross-Output = T

#how can i calculate omega with estprod?

mod2 = prodest::prodestLP(data$Y, # prodest does not have the parameter Gross-Output
             fX = data$Labor,
                 sX = data$Capital, 
             pX = data$Materials,
                 idvar = data$id,
             timevar = data$year,
                 opt='optim',
                 exit = F,
                 tol = 1e-100)  

omega = prodest::omega(mod2) #with 'prodest' I can calculate omega

Answer 1

Tras el prodtest manualmente definen omega como los residuos del modelo, es decir

$$\omega_{it} = y_{it} − (\alpha + w_{it} \beta + k_{it}\gamma)=y_{it}-\hat{y}_{it}$$

estprod no tiene ninguna documentación comprensible, pero una forma de obtener lo anterior sería extraer los coeficientes del modelo y luego calcular $y-\hat{y}$ manualmente.

Podemos hacerlo con el paquete de la escoba:

install.packages("broom")
require(broom) 
coefs<-tidy(mod1)
# A tibble: 3 x 3
  statistic   bias std.error
      <dbl>  <dbl>     <dbl>
1     -3.71   4.77      2.35
2     10.5  -11.8       4.26
3      1.73  -2.14      1.64

los nombres están apagados pero la primera columna da los mismos coeficientes que cuando inspecciono el modelo con el resumen así que son correctos.

A continuación podemos calcular $y-\hat{y}$ como

omega_1 = (data$Y - data$Labor*coefs$statistic[1] - data$Capital*coefs$statistic[2] - data$Materials*coefs$statistic[3])  

el código es poco elegante, pero debería hacer el trabajo.

Como nota al margen, podría valer la pena considerar la estimación de la PTF sin estos paquetes. Normalmente los paquetes son la mayor fortaleza de R pero en este caso estprod no tiene una buena documentación o manual que lo acompañe y el prodest no tiene la posibilidad de añadir la salida bruta todavía como en el enlace que los devs mencionaron. En casos como este puede ser mejor programar tu propia función.