Elasticity of substitution by regression: Biased results (simulation)

Question

I have the following simulation problem:

Consumers, whose utility I know, go shopping for two goods. However, prices differ each time they visit a shop. Therefore, these consumers always purchase accordingly to their relevant current budget constraints and utility functions. This gives me data for purchased quantities and corresponding prices.

Few assumptions:

• The utility is given by CES.
• Utility per consumer is stable over time.
• Neither perfect complements nor perfect substitutes.
• Coefficients for elasticity of substitution differ for consumers.

What I want to figure out from this data is the coefficient of elasticity of substitution. I tried to simulate for 1000 iterations, regress accordingly, but the results are wrong/biased.

My approach:

The problem is the following:

$$ \max_{\boldsymbol{x}} U_{i}(\boldsymbol{x}) = \left( \alpha_1 x_1^{\rho_{i}} + \alpha_2 x_2^{\rho_{i}} \right)^{\frac{1}{\rho_{i}}} \qquad s.t. \qquad P_1 x_1 + P_2 x_2 \leq M $$

Where we expect the $\rho_i$ to differ for each consumer.

I know how to compute elasticity of substitution:

$$\sigma_{ij} = \frac{\frac{\partial (x_j/x_i)}{x_j/x_i}}{\frac{\partial MRS_{ij}}{MRS_{ij}}}$$

Where $MRS_{ij} = MU_i/MU_j$.

Then, I know the main optimization condition from Lagrangian:

$$\frac{MU_i}{MU_j} = \frac{P_i}{P_j}$$

Therefore, I can plug this into the elasticity of substitution, obtaining:

$$ \sigma_{ij} = \frac{\frac{\partial (x_j/x_i)}{x_j/x_i}}{\frac{\partial P_i/P_j}{P_i/P_j}} $$

This can be rewritten in logarithmic and regression terms as:

$$\ln \left( \frac{x_j}{x_i} \right) = \sigma_{ij} \ln \left( \frac{P_i}{P_j} \right) + \epsilon$$

I know the analytic solution for this optimization problem (see code below), so I can compute the optimal purchased quantity for each consumer. I simulate prices, assuming both consumers go into the shop each day (during the day, prices are the same). Then, consumers perform their choice... I observe these choices and apply the regression equation above. However, the results are weird...

I was expecting the $\hat{\sigma_{i,j}}$ to be close to $1/(1-mean(\rho))$, however, in reality, the bigger the difference between $\rho_1$ and $\rho_2$ the bigger and more positive is the bias. It is essential to add that it gives precise results, when both parameters $\rho_1 = \rho_2$.

Now the question is: Why are those results biased? Why is the bias positive the bigger is the difference between both real elasticities? Why do I not get some mean value of elasticity of substitution? And how to properly estimate it?

The code in R used for simulation:

CES_2x_analythic = function(a1, a2, rho, P1, P2, M) {
  sc_factor = ((a2/a1)*(P1/P2))^(1/(rho-1))
  
  x1 = M/(P1+P2*(1/sc_factor))
  x2 = M/(P1*sc_factor + P2)
  
  return(c(x1, x2))
}

multicust_2goods = function(a1 = 1, a2 = 2, rho,
                            M = 100, n_iter = 1000) {
  
  n_cust = length(rho)
  
  P1 = runif(n_iter, 1, 4)
  P1 = P1[rep(1:n_iter, rep(n_cust, n_iter))]
  
  P2 = runif(n_iter, 1, 4)
  P2 = P2[rep(1:n_iter, rep(n_cust, n_iter))]
  
  rho = rho[rep(1:n_cust, n_iter)]
  mean_rho = mean(rho)
  mean_eos = 1/(1-mean_rho)
  print(mean_eos)
  
  vys = matrix(NA, nrow = n_iter*n_cust, ncol = 2)
  vys_teor = matrix(NA, nrow = n_iter*n_cust, ncol = 2)
  
  for (i in 1:n_iter) {
    
    vys[i,] = CES_2x_analythic(a1 = a1, a2 = a2, rho = rho[i],
                               P1 = P1[i], P2 = P2[i], M = M)
    
    vys_teor[i,] = CES_2x_analythic(a1 = a1, a2 = a2, rho = mean_rho,
                               P1 = P1[i], P2 = P2[i], M = M)
    
  }
  
  ln_goods = log(vys[,2]/vys[,1])
  ln_prices = log(P1/P2)
  
  OLS = lm(ln_goods ~ ln_prices -1)
  
  ln_g = log(vys_teor[,2]/vys_teor[,1])
  ln_p = log(P1/P2)
  
  OLS_teor = lm(ln_g ~ ln_p -1)
  
  return(rbind(coef(OLS), coef(OLS_teor)) )
  
}

multicust_2goods(rho = c(-0.5,+0.5), n_iter = 1000)

Thank you very much!

Answer 1

I am answering my question because I have found the solution.

The thing is: I have started with the wrong premise and therefore reached poor conclusion. Due to this I will edit the question a little bit after some time.

The problem is esentially the one of "averaging over arguments of nonlinear functions". I stated that I expected the result to be close to the $1/(1- mean(\boldsymbol\rho))$. However, this expectation is wrong and here is why:

The function of elasticity of substitution $\sigma_{i,j}$ as depending on $\rho$ is the following:

$$\sigma_{i,j} = \frac{1}{1-\rho} $$

Which might be shown in the following graph:

If we have coefficients $\rho_1$ and $\rho_2$, we can plot them in a graph (black dots), while we know what elasticity of substitution would correspond to each one of them (green line and black dots on it). Then, we can make an average of both points... BUT, here we can see, that elasticity of substitution corresponding to average of $\rho$ parameters (violet point) would always be below the value of average elasticity of substitution (orange point).

This also explains why the (pseudo-)bias I thought to see in simulation would always be positive!

When I use the correct approach, we can see that there is actually no bias at all...

So, what should the $\hat{\sigma_{i,j}}$ be equal to is this:

$$\hat{\sigma_{i,j}} = mean \left( \frac{1}{1-\rho_i} \right) $$