# Pooled model SUR

#1

I am trying to estimate a pooled SUR in R.
I am not able to replicate the results given in text.

I am also attaching my code and the example in book.

![image2|350x500]

rm(list=ls())

costdata= read.csv("TableF10-2.csv")
install.packages("micEcon")
library(micEcon)
library(tidyverse)
costdata= costdata %>%
mutate(lnpkpm= log(costdata$Pk)/log(costdata$Pm)) %>%
mutate(lnplpm= log(costdata$Pl)/log(costdata$Pm)) %>%
mutate(lnpepm= log(costdata$Pe)/log(costdata$Pm)) %>%
mutate(sk= (costdata$K*costdata$Pk)/costdata$Cost) %>% mutate(sl= (costdata$L*costdata$Pl)/costdata$Cost) %>%
mutate(se= (costdata$E*costdata$Pe)/costdata$Cost) View(costdata) library(systemfit) eqn_cpaital= sk~ lnpkpm+ lnplpm+ lnpepm eqn_labour= sl~ lnpkpm+ lnplpm+ lnpepm eqn_energy= se~ lnpkpm+ lnplpm+ lnpepm system1=list(capital= eqn_cpaital, labour= eqn_labour, energy= eqn_energy) #SUR costdata_SUR= systemfit(system1, method = "SUR", data = costdata, methodResidCov= "noDfCor", pooled = TRUE, residCovWeighted = TRUE ) summary(costdata_SUR) Any help is appreciated. Thanks! #2 #3 #4 #5 What output are you getting when you run summary(costdata_SUR)? Is the script running all the way through? It’s a bit hard to read the text from the images — but there is a nice little SUR FAQ section from UCLA’s IDRE #6 It’s running through. systemfit results method: SUR N DF SSR detRCov OLS-R2 McElroy-R2 system 72 60 0 0 0.803285 0.893478 N DF SSR MSE RMSE R2 Adj R2 capital 24 20 0 0 3.0e-05 0.629199 0.573579 labour 24 20 0 0 6.7e-05 0.794454 0.763623 energy 24 20 0 0 2.1e-05 0.920042 0.908048 The covariance matrix of the residuals used for estimation capital labour energy capital 7.69646e-10 1.19362e-09 4.26354e-10 labour 1.19362e-09 3.78213e-09 9.95445e-10 energy 4.26354e-10 9.95445e-10 3.58730e-10 The covariance matrix of the residuals capital labour energy capital 7.69646e-10 1.19362e-09 4.26354e-10 labour 1.19362e-09 3.78213e-09 9.95445e-10 energy 4.26354e-10 9.95445e-10 3.58730e-10 The correlations of the residuals capital labour energy capital 1.000000 0.699601 0.811410 labour 0.699601 1.000000 0.854604 energy 0.811410 0.854604 1.000000 SUR estimates for 'capital' (equation 1) Model Formula: sk ~ lnpkpm + lnplpm + lnpepm Estimate Std. Error t value Pr(>|t|) (Intercept) 2.74321e-04 4.93844e-05 5.55481 1.9469e-05 *** lnpkpm 1.31729e-05 5.60644e-06 2.34960 0.0291753 * lnplpm -8.49374e-05 2.60368e-05 -3.26220 0.0039001 ** lnpepm 5.62066e-05 8.91539e-06 6.30445 3.7346e-06 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 3e-05 on 20 degrees of freedom Number of observations: 24 Degrees of Freedom: 20 SSR: 0 MSE: 0 Root MSE: 3e-05 Multiple R-Squared: 0.629199 Adjusted R-Squared: 0.573579 SUR estimates for 'labour' (equation 2) Model Formula: sl ~ lnpkpm + lnplpm + lnpepm Estimate Std. Error t value Pr(>|t|) (Intercept) 1.46879e-03 1.09474e-04 13.41675 1.8472e-11 *** lnpkpm -1.76772e-05 1.24282e-05 -1.42234 0.170337 lnplpm -1.59515e-04 5.77179e-05 -2.76370 0.011979 * lnpepm 1.45536e-04 1.97634e-05 7.36390 4.0910e-07 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 6.7e-05 on 20 degrees of freedom Number of observations: 24 Degrees of Freedom: 20 SSR: 0 MSE: 0 Root MSE: 6.7e-05 Multiple R-Squared: 0.794454 Adjusted R-Squared: 0.763623 SUR estimates for 'energy' (equation 3) Model Formula: se ~ lnpkpm + lnplpm + lnpepm Estimate Std. Error t value Pr(>|t|) (Intercept) 3.52891e-04 3.37154e-05 10.46675 1.4612e-09 *** lnpkpm -1.37385e-05 3.82759e-06 -3.58932 0.0018333 ** lnplpm -1.29858e-04 1.77757e-05 -7.30536 4.6046e-07 *** lnpepm 7.39568e-05 6.08666e-06 12.15065 1.0902e-10 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 2.1e-05 on 20 degrees of freedom Number of observations: 24 Degrees of Freedom: 20 SSR: 0 MSE: 0 Root MSE: 2.1e-05 Multiple R-Squared: 0.920042 Adjusted R-Squared: 0.908048 #7 So it’s the actual values that are coming out wrong? Sorry— I honestly don’t know much about SUR, so I couldn’t tell if you were having a code issue or a calculation one! #8 Yes the actual values are wrong. #9 It’s been a long time since I read Greene, and only ever ran one SUR (in Stata, not R), so I could be completely off. First, I think the log quantities in your code should be lnpkpm = log(Pk / Pm) instead of lnpkpm = log(Pk) / log(Pm) given the equations. Next, I think K, L, and E are already factor shares, according to the descripton of those variables on the dataset page. Finally, I think the SUR estimates are in Table 10.3; Table 10.4 involves using those estimates and fitted cost shares to get the relevant elasticities. That said, I tried my corrections and got estimates similar in magnitude and direction to those in Table 10.3, but not exactly those numbers. I’m not sure why…it might mean I’m totally off, so I would take what I say here with a pinch of salt, and defer to those with more experience with the translog method. #10 Apologies for not posting table 10.3. The estimates are still not matching, even after emploing the suggestions. Any further help is appreciated. Thanks. #11 Nobody has any suggestion?? #12 You’re outside my domain of expertise! I’m not sure if the economics Stack Exchange site is “mature” yet, but you might give this question a go either there, on Cross-Validated: There’s also a Principles of Econometrics in R bookdown by Constantin Colonescu, that has an SUR model in this example. If you can break it down to a specific piece of code that isn’t evaluating as expected, you might have a better chance getting help here! Hopefully these suggestions can meta-help you find the right help! #13 Thanks, Happy Birthday!! #14 Howdy. Other than the setup, this isn’t really an R or RStudio question. But happy to help. It’s been 8 years for me spending much time with Greene The model setup by eq 10-40 is, I think: costdata = read.csv(".../TableF10-2.csv") # setup dependent vars for /delta_... costdata = costdata %>% mutate( lnpkpm = log(Pk / Pm), lnplpm = log(Pl / Pm), lnpepm = log(Pe / Pm) ) setup eq 10-40: df_pooledreg <- bind_rows( tibble( y = costdata$K,
k = 1,
l = 0,
e = 0,
kk = costdata$lnpkpm, kl = costdata$lnplpm,
ke = costdata$lnpepm, ll = 0, le = 0, ee = 0, year = costdata$Year
),
tibble(
y = costdata$L, k = 0, l = 1, e = 0, kk = 0, kl = costdata$lnpkpm,
ke = 0,
ll = costdata$lnplpm, le = costdata$lnpepm,
ee = 0,
year = costdata$Year ), tibble( y = costdata$E,
k = 0,
l = 0,
e = 1,
kk = 0,
kl = 0,
ke = costdata$lnpkpm, ll = 0, le = costdata$lnplpm,
ee = costdata$lnpepm, year = costdata$Year
)
)

run “The Pooled Model” eq 10-19 in Greene chapter 10.
Note, eq 10-19 doesn’t have an intercept.

pooledreg = glm(
y ~ 0 + k + l + e + kk + kl + ke + ll +le + ee,
data = df_pooledreg
)
summary(pooledreg)

I get the following coefficients:

Coefficients:
Estimate Std. Error t value Pr(>|t|)
k   0.056223   0.001888  29.773  < 2e-16 ***
l   0.253409   0.001852 136.834  < 2e-16 ***
e   0.041861   0.002291  18.272  < 2e-16 ***
kk  0.030426   0.008067   3.772 0.000349 ***
kl  0.001755   0.004556   0.385 0.701342
ke -0.003562   0.007481  -0.476 0.635562
ll  0.075127   0.005280  14.228  < 2e-16 ***
le  0.003252   0.005756   0.565 0.574043
ee  0.046777   0.017415   2.686 0.009136 **
---

Now, Table 10.4. The first line, Fitted shares is just:

fitted_k = predict(
pooledreg,
df_pooledreg %>% filter(year == 1959 & k)
)
fitted_l = predict(
pooledreg,
df_pooledreg %>% filter(year == 1959 & l)
)
fitted_e = predict(
pooledreg,
df_pooledreg %>% filter(year == 1959 & e)
)
# Note that assuming a constant returns to scale production function,
# m is just 1-(k+l+e)
fitted_m = 1 - (fitted_k + fitted_l + fitted_e)

fitted_k: 0.05655583
fitted_l: 0.2747627
fitted_e: 0.04487152
fitted_m: 0.6238099

How do we estimate the derivatives involving Materials? i.e. those $/delta_{M...}$? For that, turn to equation 10-38 (well, it’s eq 10-38 in my Greene, 6th edition of econometric analysis anyway. they seem to be moving things around). That equation just says,
$\Sigma_{i=1}^M \delta_{ij} = 0$ and $\Sigma_{j=1}^M \delta_{ij} = 0$.

Really simple, you can see this in table 10.4. It just says that $\delta_{KK} +$\delta_{KL} + $\delta_{KE} +$\delta_{KM} = 0.
You’ll notice in table 10.4 any derivative with an M has an asterisks next to it, referring to eq 10-38.

Equation 10-39 gives the own and cross elasticities of substitution.
So with $\theta_{kk}$ for example, (/delta_{kk} + s_{i}(s_{i} - 1)) / s_{i}^2, where i is capital share’s output in that year.
or with my estimates and your data:

s_k1959 = ((costdata %>% filter(Year == 1959)))$K delta_kk = pooledreg$coefficients['kk'] %>% as.vector
(delta_kk + s_k1959 * (s_k1959 - 1)) / (s_k1959^2)
[1] -7.214395 # that's under "Implied Elasticities of Substitution, 1959" for Capital/Capital

And for cross price elasticity of substitution, eq 10-39 says
$\theta_{ij}$ is , (/delta_{ij} + s_{i} * s{j}) / (s_{i} * s_{j}),
for example for $\theta_{KL}$:

s_k1959 = ((costdata %>% filter(Year == 1959)))$K s_l1959 = ((costdata %>% filter(Year == 1959)))$L
delta_kl = pooledreg\$coefficients['kl'] %>% as.vector
(delta_kl + (s_k1959 * s_l1959)) / (s_k1959 * s_l1959)
[1] 1.103912 # this is our estimate for "Implied Elasticities of Substitution, 1959" for Labor/Capital

There are slight deviations. They might come from using an updated dataset…?

#15

Thanks Curtis the values are indeed slightly different with my working also. I just want to know why the pooled SUR from systemfit fails here??

#16

Does the systemfit package throw an error/the code itself fail somehow? If so, you might consider reporting it as a bug or feature request on systemfit's R-Forge page. Here’s the direct link to its trackers.

#17

No the output produced is different and the code runs through. I am already intouch with the developers there, Thanks.