Preliminary Analysis: Are franchisees controlled or infleuenced by franchisors?

A structual approach is used to test the following hypothesis: Franchisors have high levels of control over franchisees in metropolitan cities?

I empirically estimate the level of control of franchisors over franchisees by using the demand and supply estimation.

Demand Estimation

The indirect utility function that consumers who purchase product \(j\) in market \(t\) is

\[\begin{array} {lll} u_{ijt} & = & \alpha p_{jt} + X_{jt}\beta + \xi_{j} + \zeta_{ig} + (1-\sigma)\epsilon_{ijt} \\ & = & \delta_{jt} + \zeta_{ig} + (1-\sigma)\epsilon_{ijt} \end{array} \] - \(g=0, 1,\cdots, G\) represents a group in which a product is associated. - \(F_g\) is a set of products in group \(g\). - \(\zeta\) is common for all products in group \(g\) - \(\epsilon\) is IID with the extreme value distribution.

With a dummy variable \(d_{jg}\) (it is one if $j F_g $), the utility function can be rewritten:

\[ u_{ijt} = \delta_{jt} + \sum_g[d_{jg} \zeta_{ig}] + (1-\sigma)\epsilon_{ijt} \]

The market share of product \(j\) (the probability that consumers choose product \(j\)) is

\[ s_{jt} = s_{jgt} \cdot s_{g} \] where \(s_{jg}\) represents the probability that product \(j\) is chosen given all products in \(g\), while \(s_g\) is the probability that group \(g\) is chosen.

\(s_jg\) is defined as follows: \[ s_{jg} = \frac{\exp(\frac{\delta_j}{1-\sigma})}{\sum_{j\in F_g}\exp(\frac{\delta_j}{1-\sigma})} = \exp(\frac{\delta_j}{1-\sigma}) / D_g \] where \[ D_g = \sum_{j\in F_g} \exp(\frac{\delta_j}{1-\sigma}) \]

\(s_g\) is the following:

\[ s_g = \frac{D_g^{(1-\sigma)}}{\sum_{g'}D_{g'}^{(1-\sigma)}} \] Market share of product \(j\) is

\[ s_{j} = s_{jg} \cdot s_g = \frac{\exp(\frac{\delta_j}{1-\sigma})}{D_g^{\sigma}[\sum_{g'}D_{g'}^{(1-\sigma)}]} \]

Outside option is defined as follows: \[ s_0 = \frac{1}{\sum_{g'}D_{g'}^{(1-\sigma)}} \]

Take the log of market share of product \(j\) and within group market share of product \(j\),

\[ \begin{align} \ln(s_j) - \ln(s_0) & = \delta_j / (1-\sigma) - \sigma \ln(D_g) && Eq. 1 \\ \ln(s_{jg}) & = \delta_j / (1-\sigma) - \ln(D_g) && Eq. 2 \end{align} \]

Substracting \(\sigma \cdot Eq. 2\) from \(Eq. 1\), the equation is rewritten as follows:
\[ \ln(s_{jt}) - \ln(s_{0t}) = \alpha p_{jt} + x_{jt} \beta + \sigma\ln(s_{jg}) + \epsilon_{jt} \]

Price Elasticities and \(\partial s / \partial p\)

There are three different formulat for price elasticities under the nested logit model: $$ \[\begin{eqnarray} \eta_{jk} & = & \frac{\partial s_{j}}{\partial p_{k}} \frac{p_{k}}{s{j}} \\ & = & [s_j \cdot \frac{\alpha}{1-\sigma}(1- \sigma s_{jg} -(1-\sigma)s_j )] \frac{p_j}{s_j} && (j = k) \\ & = &[- s_j \cdot\frac{\alpha}{1-\sigma} (\sigma s_{kg} + (1-\sigma)s_k )] \frac{p_k}{s_j} && (j \neq k; j,k \in g)] \\ & = & [-s_j \cdot \alpha \cdot s_k ] \frac{p_k}{s_j} && (j\neq k; j\in g, k\notin g ) \end{eqnarray}\] $$

When calculating \(\partial s / \partial p\), the estimated market shares (\(s_j, s_{jg}\)) are used, and $<0 $.

Market Definition

Markets are defined based on the cluster analysis. I use the density based clustering analysis with noise (DBSAN). For detail of analysis, see http://rpubs.com/jhkoh17/488032.

Supply Side and Estimation of the Level of Control from Franchisors

The profit function of a firm is as follows:

\[\begin{eqnarray} \Pi_{jt} & = & \underbrace{(p_{jt} - mc_{jt})s_{jt}M_{t}}_{\text{own profits}(\pi_{jt})} + \underbrace{\sum_{k\neq j}\lambda(p_{kt}-mc_{kt})s_{kt}M_t}_{\text{considering others' profits}} \end{eqnarray}\]

where \(\lambda\) is a \((J-1) \times (J-1)\) matrix measuring the level of controls. This would vary depending on assumptions of market comeptition and/or vertical relationship between franchisors and franchisees.

If firm \(j\) only considers, \(\lambda\) would be a zero matrix and \(\Pi_{jt} = \pi_{jt}\).

However, firm \(j\) considers and is influenced by other firms under the same franchisors, \(\lambda\) is non-zero matrix. This equation can be rewritten in a matrix form:

\[\begin{eqnarray} \Pi = \Lambda(P - MC)SM_t \end{eqnarray}\]

where \(\Lambda\) is a \(J \times J\) matrix with ones in its diagonal elements.

Studies of mergers and product differentiation assumes the

Case 1: Franchisor influence

For example, there are three firms in the markets (Firm 1, 2, and 3). Firms 1 and 2 are under the same franchisor, while Firm 3 is associated with another franchisor. Thus,

\[\Lambda^1 = \begin{bmatrix} 1 & \lambda^1 & 0\\ \lambda^1 & 1 & 0 \\ 0 & 0 & 1 \\ \end{bmatrix} \]

Case 2: Frnachisors and collusion with other competitiors

This matrix can be expended if collusion among firms occurs as follows:

\[\Lambda^{12} = \begin{bmatrix} 1 & \lambda^1 & \lambda^2\\ \lambda^1 & 1 & \lambda^2 \\ \lambda^2 & \lambda^2 & 1 \\ \end{bmatrix} \] I assume \(\lambda^1\) is the same for all firms under the same franchisors, while \(\lambda^2\) is the same for all competitiors.

Here I focus on Case 1 and move onto Case 2.

Model for Marginal Costs

The model for marginal costs is used to create the moment conditions for GMM for the supply side, like BLP (1995). The model is as follows:

\[ mc_{jt} = w_{jt}\gamma + \omega_{jt} \] where \(mc\) is derived from the demand side model, \(w\) is a matrix of factors affecting costs.

The first order condition of Case 1

\[\begin{eqnarray} \frac{\partial \Pi_jt}{\partial p_jt} & = & 0 \end{eqnarray}\] \[ s_{jt} + (p_{jt}-mc_{jt})\frac{\partial s_jt}{\partial p_{jt}} + [\sum_{k} \lambda_{jk} (p_{kt} - mc_{kt})\frac{\partial s_{kt}}{ \partial p_{jt}}] = 0 \] \[ s_{jt} + (p_{jt}-mc_{jt})\frac{\partial s_{jt}}{\partial p_{jt}} + [\sum_{k} \lambda^1 (p_{kt} - mc_{kt})\frac{\partial s_{kt}}{ \partial p_{jt}}] = 0 \] where \(j\) and \(k\) are under the same franchisor.

Define \(\Omega = -\Lambda^1 \cdot \partial s / \partial p\)( \(\partial s / \partial p\) is a matrix form of \(\frac{\partial s_{jt}}{\partial p_{kt}}\)) and rewrite the above equation is the marginal cost fucntion

\[ mc(\lambda) = p - \Omega^{-1}(\lambda)\cdot s() \] With the function of \(mc\), the shock for the marginal cost function:

\[ \omega = p- w\gamma - \Omega^{-1}()\cdot s() \]

With the instruments for the supply side, \(Z_s\), the moment condits are defined as \[ g(\theta)=E(Z_s\omega) = 0 \] where \(\theta = [\lambda, \gamma]\).

Using the moment conditions, I estimate two set of parameters, \(\gamma\) and \(\lambda\). \(\lambda \in(0,1)\) is similar to conduct parameters used in the literature in multi-market contact (See Ciliberto and Williams(2014))

Estimation Detail

The estimated marginal costs from the supply side is the following: \[ \hat{mc} = W \cdot \hat{\gamma} \]

Rewrite the above the equation: \[\begin{align} \omega & = mc - \hat{mc} \\ & = p - \Omega^{-1}() s() - \hat{mc} \\ \end{align}\]

The sample momoent conditions are as follows: \[\begin{align} \hat{g(\theta)} & = \sum_i^n Z_i (mc_i - W_i\hat{\gamma})/n \\ & = Z' (mc - W\gamma) / n \end{align}\] Finally, the GMM objective function is the following: \[ \arg \min_\theta \hat{g}(\theta)' A \hat{g}(\theta) \] where \(A\) is a weighting matrix.

Demand Estimation

The nested logit demand estimation

# Nested Logit Model ------------------------------------------------------
# y.fld 
load(file="df_hou4.Rdata")
load(file="blp_mmc.Rdata")
y.fld <- "y2"
# # prc.fld 
end.fld <- "adr"
# x.fld
x.fld <- c( "n_todo", "n_room_amenity", "n_service",  "cbd",  "air", "log(sjg2)")

# iv.fld
iv.fld <- colnames(blp)
iv.fld <- iv.fld[c(3,4,7,8,11,12, 13,14, 15, 16)]

# nest.var.fld
nest.var.fld <- c("log(sjg2)")

ols.eq <- paste0(y.fld, " ~ ", end.fld, " + ", paste(x.fld, collapse = " + "))
ols.eq
## [1] "y2 ~ adr + n_todo + n_room_amenity + n_service + cbd + air + log(sjg2)"
iv.eq <- paste0(ols.eq, " | ", paste(x.fld, collapse = " + "), " + ", paste(iv.fld, collapse=" + "))
iv.eq
## [1] "y2 ~ adr + n_todo + n_room_amenity + n_service + cbd + air + log(sjg2) | n_todo + n_room_amenity + n_service + cbd + air + log(sjg2) + sum_same_n_todo + sum_other_n_todo + sum_same_n_service + sum_other_n_service + sum_same_n_room_amenity + sum_other_n_room_amenity + sum_same_cbd + sum_other_cbd + sum_same_air + sum_other_air"
nl.est <- ivreg(iv.eq, data=df_hou4)
summary(nl.est)
## 
## Call:
## ivreg(formula = iv.eq, data = df_hou4)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -4.08924 -0.47677  0.04444  0.56038  4.51659 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     0.649766   0.136789   4.750 2.23e-06 ***
## adr            -0.015212   0.001501 -10.134  < 2e-16 ***
## n_todo          0.114472   0.023562   4.858 1.31e-06 ***
## n_room_amenity -0.032589   0.020432  -1.595 0.110920    
## n_service       0.060663   0.017722   3.423 0.000636 ***
## cbd             1.194139   0.182347   6.549 7.93e-11 ***
## air            -0.087959   0.140836  -0.625 0.532360    
## log(sjg2)       0.919510   0.024590  37.394  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.062 on 1513 degrees of freedom
## Multiple R-Squared: 0.3366,  Adjusted R-squared: 0.3335 
## Wald test: 228.5 on 7 and 1513 DF,  p-value: < 2.2e-16

Obtaining Market Share

sj.nl <- sj.fn(est=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)", 
               nest.fld="class2", mkt.id.fld="mkt.id",dat = df_hou4)

Jacobian of Market Share with the respect to price

# dsdp --------------------------------------------------------------------
dshare.dp <- function(est.data, prc.fld, nest.var.fld, nest.fld, mkt.id.fld, brand.id.fld, dat, lambda, ...){
  # data prepration 
  est <- est.data$coefficient
  xi <- est.data$residuals
  mkt.id <- as.matrix(dat[ , mkt.id.fld])
  brand.id <- as.matrix(dat[, brand.id.fld])
  nest <- as.matrix(dat[, nest.fld])
  beta <- as.matrix(est[!names(est) %in% nest.var.fld])
  alpha  <- beta[prc.fld,]
  sigma  <- est[nest.var.fld]
  sj.dat <- sj.fn(est = est.data, prc.fld=prc.fld, nest.var.fld=nest.var.fld, 
                  nest.fld=nest.fld, mkt.id.fld=mkt.id.fld, dat=dat)
  print("share calculated")
  
  # sj  <- as.matrix(sj.dat[,1])
  # sjg <- as.matrix(sj.dat[,2])
  
  sj  <- as.matrix(dat[, "sj2"])
  sjg <- as.matrix(dat[, "sjg2"])
  
  
  sj.ob <- as.matrix(dat[, "sj2"])
  
  
  # three cases 
  # case 1: own price elasticities 
  # case 2: cross price elastiticies within nest
  # case 3: cross price elasticities outside nest (the same as in logit model)
  dsdp.list <- list()
  dsdp.ow.list <- list()
  
  for (t in unique(mkt.id)){
    inc <- mkt.id == t
    nest.mkt <- nest[inc]
    sj.mkt <- sj[inc]
    sjg.mkt <- sjg[inc]
    sj.ob.mkt <- sj.ob[inc]
    nest.mkt.mat <- outer(nest.mkt, nest.mkt, function(x,y) (x ==y)*1) # matrix of nest 
    part1 <- matrix(0, nrow=length(sj.mkt), ncol=length(sj.mkt), byrow=T)
    diag(part1) <- 1/(1-sigma) * sj.mkt
    part2 <- sigma /(1-sigma) * tcrossprod(sj.mkt, sjg.mkt)
    part3 <- tcrossprod(sj.mkt, sj.ob.mkt)
    dsdp.mkt <- alpha*(part1 - part2 - part3)
    brand.id.mkt <- brand.id[inc]
    owner.mat <- outer(brand.id.mkt, brand.id.mkt, function(x,y) x==y) * lambda
    diag(owner.mat) <- 1
    dsdp.ow.mkt <- alpha*(part1 - part2 - part3) * owner.mat
    # dsdp.mkt <- - alpha * tcrossprod(sj.mkt, sj.mkt) # out
    # dsdp.mkt <- (-(alpha * sigma /(1-sigma)) * tcrossprod(sj.mkt, sjg.mkt)) * nest.mkt.mat + dsdp.mkt # cross in the same nest
    # dsdp.mkt <- diag(alpha * sigma /(1-sigma) * sj.mkt) + dsdp.mkt # own 
    
    dsdp.list[[t]] <- dsdp.mkt
    dsdp.ow.list[[t]] <- dsdp.ow.mkt
  }
  return(dsdp.ow.list)
}
# names(data.houston.cl.iv1)

dsdp1 <- dshare.dp(est.data=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)", 
                   nest.fld="class2", mkt.id.fld="mkt.id", brand.id.fld="id.brand2",
                   lambda = 1, dat = df_hou4)
## [1] "share calculated"
dsdp2 <- dshare.dp(est.data=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)", 
                  nest.fld="class2", mkt.id.fld="mkt.id", brand.id.fld="id.brand2", 
                  lambda = 0.8, dat = df_hou4)
## [1] "share calculated"
# Marginal costs ----------------------------------------------------------
mc <- function(est.data, prc.fld, nest.var.fld, nest.fld, mkt.id.fld, brand.id.fld, lambda = 1, dat, ...){
  # data prepration 
  
  est <- est.data$coefficient
  xi <- est.data$residuals
  mkt.id <- as.matrix(dat[ , mkt.id.fld])
  brand.id <- as.matrix(dat[, brand.id.fld])
  # 
  brand.id[brand.id==1] <- (max(brand.id)+1):(max(brand.id)+sum(brand.id==1))
  price <- as.matrix(dat[, prc.fld])
  nest <- as.matrix(dat[, nest.fld])
  beta <- as.matrix(est[!names(est) %in% nest.var.fld])
  alpha  <- beta[prc.fld,]
  sigma  <- est[nest.var.fld]
  # dsdp <- dshare.dp(est=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)", 
  #                    nest.fld="class2", mkt.id.fld="mkt.id", brand.id.fld="id.brand2", 
  #                    lambda = 1, dat = df_hou4)
  # 
  
  dsdp <- dshare.dp(est.data =  est.data, prc.fld=prc.fld,
                    nest.var.fld=nest.var.fld, 
                    nest.fld=nest.fld, mkt.id.fld=mkt.id.fld,
                    brand.id.fld=brand.id.fld, 
                    lambda = lambda, dat = dat)
  
  
  print("dshare.dprice calcuated")
  
  sj.dat <- sj.fn(est = est.data, prc.fld=prc.fld, nest.var.fld=nest.var.fld, 
                  nest.fld=nest.fld, mkt.id.fld=mkt.id.fld, dat=dat)
  print("share calculated")
  sj  <- as.matrix(sj.dat[,1])
  sjg <- as.matrix(sj.dat[,2])
  
  # sj <- dat[, "sj"]
  # sjg <- dat[, "sjg"]
   
  mc.list <- list()
  for(t in unique(mkt.id)){
    inc <- mkt.id == t
    price.mkt <- price[inc]
    sj.mkt <- sj[inc]
    dsdp.mkt <- dsdp[[t]] * -1
    brand.id.mkt <- brand.id[inc]
    mc.mkt <- price.mkt - solve(dsdp.mkt) %*% sj.mkt
    mc.list[[t]] <- mc.mkt
  }
  return(unlist(mc.list))
}
mc1 <- mc(est.data=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)",
          nest.fld="class2", mkt.id.fld="mkt.id", brand.id.fld="id.brand2", 
          lambda = 0.5, dat=df_hou4)
## [1] "share calculated"
## [1] "dshare.dprice calcuated"
## [1] "share calculated"
sum(mc1<0)
## [1] 107
mc2 <- mc(est.data=nl.est, prc.fld="adr", nest.var.fld="log(sjg2)", nest.fld="class2", 
          mkt.id.fld="mkt.id", brand.id.fld="id.brand2", lambda = 1, dat=df_hou4)
## [1] "share calculated"
## [1] "dshare.dprice calcuated"
## [1] "share calculated"
sum(mc2<0)
## [1] 108

Supply Side

chain.dummy.fld <- names(df_hou4)[grep("chain", names(df_hou4))][7:30]

w.fld <- c("const", "room", "s_rating", "n_todo", "n_room_amenity", chain.dummy.fld)
w.fld <- c("const", "room", "class2", "n_todo", "n_room_amenity", "factor(brand)")
w.fld <- c("const", "room", "n_todo", "n_room_type", "n_room_amenity", "n_service", "factor(brand)")
w.fld <- c("const", "room", "n_todo", "n_room_type", "n_room_amenity", "n_service", "factor(chain)")
w.fld <- c("const", "room", "n_todo", "n_room_type", "n_room_amenity", "n_service", "factor(chain)")
w.fld <- c("const",  "n_todo", "n_room_amenity", "n_service",  "cbd",  "air", "n.same.brand", "n.same.chain")


w.fld <- c("const", "room", "n_todo", "n_room_amenity", "n_service",  "cbd",  "air", "factor(chain)")
w.fld <- c("const",  "s_rating", "n_todo", "n_room_amenity", "n_service",  "cbd",  "air", "factor(brand)")

w.fld <- c("const", "room",  "s_rating", "n_todo", "n_room_amenity", "n_service",  "cbd",  "air", "factor(brand)")

# The below specification is the best specification
w.fld <- c("const", "room", "n_room_amenity", "n_room_type", "n_service", chain.dummy.fld)
ols.supply <- paste0("mc1 ~ ", paste(w.fld[-1], collapse = " + "))
ols.supply
## [1] "mc1 ~ room + n_room_amenity + n_room_type + n_service + chain1 + chain2 + chain3 + chain4 + chain5 + chain7 + chain11 + chain12 + chain13 + chain14 + chain15 + chain16 + chain17 + chain18 + chain19 + chain22 + chain23 + chain24 + chain26 + chain27 + chain28 + chain29 + chain30 + chain40"
est.supply <- lm(ols.supply, data =df_hou4)
summary(est.supply)
## 
## Call:
## lm(formula = ols.supply, data = df_hou4)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2203.84    -8.48    13.71    37.07   348.68 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -0.24398    7.46599  -0.033  0.97394    
## room             0.19881    0.03278   6.066 1.66e-09 ***
## n_room_amenity   1.53320    2.84700   0.539  0.59029    
## n_room_type     -0.67819    4.93005  -0.138  0.89061    
## n_service       -1.79414    2.12062  -0.846  0.39766    
## chain1         102.36169   12.02989   8.509  < 2e-16 ***
## chain2         106.26326   16.05974   6.617 5.10e-11 ***
## chain3          83.82054   13.47421   6.221 6.41e-10 ***
## chain4          26.97229   14.94558   1.805  0.07132 .  
## chain5          60.05010   13.46375   4.460 8.80e-06 ***
## chain7          10.51704   30.25752   0.348  0.72820    
## chain11         15.95243   22.51450   0.709  0.47872    
## chain12         68.67512   22.95114   2.992  0.00281 ** 
## chain13          6.25523   44.78437   0.140  0.88894    
## chain14         35.75558   45.26771   0.790  0.42973    
## chain15         59.62543   37.14829   1.605  0.10869    
## chain16        101.93693   22.11358   4.610 4.38e-06 ***
## chain17         40.16555   18.82108   2.134  0.03300 *  
## chain18        192.20576   63.67495   3.019  0.00258 ** 
## chain19        380.33937  125.12916   3.040  0.00241 ** 
## chain22          5.40600   23.39219   0.231  0.81727    
## chain23         46.70469   17.84744   2.617  0.00896 ** 
## chain24         15.87384   21.85492   0.726  0.46775    
## chain26         93.62352   45.81261   2.044  0.04117 *  
## chain27         15.74198   29.35814   0.536  0.59190    
## chain28         53.35809   65.11541   0.819  0.41267    
## chain29         -5.00702   44.74111  -0.112  0.91091    
## chain30         12.11475   37.32690   0.325  0.74556    
## chain40         -7.12241   15.17018  -0.470  0.63878    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 124.7 on 1492 degrees of freedom
## Multiple R-squared:  0.1691, Adjusted R-squared:  0.1535 
## F-statistic: 10.84 on 28 and 1492 DF,  p-value: < 2.2e-16

Moment conditions

est.data      <- nl.est
prc.fld       <- "adr"
nest.var.fld  <- "log(sjg2)"
nest.fld      <- "class2" 
mkt.id.fld    <- "mkt.id"
brand.id.fld  <- "id.brand2"
lambda <- 1
dat <- df_hou4

moments <- function(lambda, est.data, prc.fld, nest.var.fld, nest.fld, mkt.id.fld, 
                    brand.id.fld, w.fld, iv.fld, mc0, 
                    dat, ...){
  
  # GMM momment conditiosn
  # Parameters to be estimated: lambda and gamma 
  
  # data prepration 
  est <- est.data$coefficient # coefficients from demeand side
  mkt.id <- as.matrix(dat[ , mkt.id.fld]) # 
  brand.id <- as.matrix(dat[, brand.id.fld]) # brand.id
  w  <- as.matrix(ones=rep(1, nrow(dat)), dat[w.fld])   # RHS variable in MC: MC = w*gamma + omega
  iv <- as.matrix(dat[iv.fld]) # instruments 
  z <- cbind(w, iv) # Z = {w, iv}
  # 
  brand.id[brand.id==1] <- (max(brand.id)+1):(max(brand.id)+sum(brand.id==1))
  price <- as.matrix(dat[, prc.fld])
  nest <- as.matrix(dat[, nest.fld])
  beta <- as.matrix(est[!names(est) %in% nest.var.fld]) # coefficients of price and product attributes
  alpha  <- beta[prc.fld,]  # price coefficients
  sigma  <- est[nest.var.fld]   # nest coefficient
  # dsdp
  
  dsdp <- dshare.dp(est.data =  est.data, prc.fld=prc.fld,
                    nest.var.fld=nest.var.fld, 
                    nest.fld=nest.fld, mkt.id.fld=mkt.id.fld,
                    brand.id.fld=brand.id.fld, 
                    lambda = lambda, dat = dat)
    
  sj.dat <- sj.fn(est = est.data, prc.fld=prc.fld, nest.var.fld=nest.var.fld, 
                  nest.fld=nest.fld, mkt.id.fld=mkt.id.fld, dat=dat)
  print("share calculated")
  sj  <- as.matrix(sj.dat[,1])
  sjg <- as.matrix(sj.dat[,2])
  
  # initial value of gamma
  gamma0 <-  solve(t(w) %*% w) %*% t(w) %*% mc0
  theta0 <- cbind(lambda, t(gamma0))
  
  # updated weighting matrix
  omega.hat <- mc0 - w%*%gamma0
  z.hat <- z * matrix(rep(omega.hat, ncol(z)), ncol = ncol(z))
  
  gmm_obj <- function(theta){
    lambda <- theta[1]
    gamma <- theta[-1]
    mc <- mc(est.data = est.data, prc.fld = prc.fld, nest.var.fld = nest.var.fld, 
             nest.fld = nest.fld, mkt.id.fld = mkt.id.fld, brand.id.fld = brand.id.fld,
             lambda = lambda, dat = dat)
    
    # Projection matrix with solve(t(Z)'Z) (2sls)
    PZ <- z %*% solve(t(z) %*% z) %*% t(z)
    # updating weighting matrix 
    omega.hat <- mc - w %*% gamma
    z.hat <- z * matrix(rep(omega.hat, ncol(z)), ncol = ncol(z))
    a.inv <- try(solve(t(z.hat) %*% z.hat), silent =FALSE)
    if("matrix" == class(a.inv)){
      print("GMM step 2 updated")
      PZ <<- z %*% a.inv %*% t(z);
      PX.inv <- solve(t(w) %*% PZ %*% w)
      # theta1 <<- PX.inv %*% t(X) %*% PZ %*% delta
      # xi.hat <<- delta - X %*% theta1
      # X.hat <- (PZ %*% X) * matrix(rep(xi.hat, K), ncol = K)
      # tsls.se <- sqrt(diag(PX.inv %*% t(X.hat) %*% X.hat %*% PX.inv))
      # print("GMM step 2 updated theta1 estimate:")
      # print(beta.est <<- list ( linear = data.frame(beta.est = theta1, beta.se = tsls.se),nonlinear = theta2))
    }
    dat[, "omega.hat"] <<- omega.hat
    f <- t(omega.hat) %*% PZ %*% omega.hat
    print("Updated GMM objective:")
    print(f <- as.numeric(f))
return(f)
}
  # # 2sls for gamma
  # str.ivreg.y  <- "mc ~ "
  # str.ivreg.x  <- paste(w.fld, collapse = " + ")
  # str.ivreg.iv <- paste(iv.fld, collapse = " + ")
  # print("2SLS specification:")
  # print(fm.ivreg <- paste0(str.ivreg.y, str.ivreg.x , " | ", str.ivreg.x, " + ", str.ivreg.iv))
  # rm(str.ivreg.y, str.ivreg.x, str.ivreg.prc, str.ivreg.iv)
  # mo.ivreg <- ivreg(fm.ivreg, data = dat, x = TRUE))
  
  # stop condition
  tol <- 1e-8
  max.it <- 1000
  # lower and upper bounds
  lb <- c(0, rep(-Inf, ncol(z)))
  up <- c(1, rep(Inf, ncol(z)))
  
  # optimization 
  theta.est <- neldermead(x0 = theta0, fn = gmm_obj, lower=NULL, upper=NULL, 
                    control = list(xtol_rel = tol, maxeval = max.it))
  theta <- theta.est$par
  gmm.obj.value <- gmm_obj(theta)/nrow(dat)
  
return(list(theta=theta, gmm.obj.value = gmm.obj.value))
}
std.err <- function(lambda, est.data, prc.fld, nest.var.fld, nest.fld, mkt.id.fld, 
                    brand.id.fld, w.fld, iv.fld, mc0, 
                    dat, ...){
  
  
  
  # Data Prep
  est <- est.data$coefficient # coefficients from demeand side
  mkt.id <- as.matrix(dat[ , mkt.id.fld]) # 
  brand.id <- as.matrix(dat[, brand.id.fld]) # brand.id
  w  <- as.matrix(ones=rep(1, nrow(dat)), dat[w.fld])   # RHS variable in MC: MC = w*gamma + omega
  iv <- as.matrix(dat[iv.fld]) # instruments 
  z <- cbind(w, iv) # Z = {w, iv}
  # 
  brand.id[brand.id==1] <- (max(brand.id)+1):(max(brand.id)+sum(brand.id==1))
  price <- as.matrix(dat[, prc.fld])
  nest <- as.matrix(dat[, nest.fld])
  beta <- as.matrix(est[!names(est) %in% nest.var.fld]) # coefficients of price and product attributes
  alpha  <- beta[prc.fld,]  # price coefficients
  sigma  <- est[nest.var.fld]   # nest coefficient
  # dsdp
  
  
  # dshare.dp with lambda == 1
    dsdp <- dshare.dp(est.data =  est.data, prc.fld=prc.fld,
                    nest.var.fld=nest.var.fld, 
                    nest.fld=nest.fld, mkt.id.fld=mkt.id.fld,
                    brand.id.fld=brand.id.fld, 
                    lambda = 1, dat = dat)

      
}

Estimation Result

  • \(\lambda\) : Estimated Conduct Parameter
  • Others : Coefficients of the Supply Side \(mc = W\gamma + \omega\)
##               var          coefficient
## 1          Lambda    0.305436988717096
## 2        Constant     82.6197969408459
## 3            Room -0.00132904231676757
## 4  n_room_amenity    -1.55471989587518
## 5     n_room_type    0.249417985125057
## 6       n_service     1.35556165353253
## 7          chain1    -7.11104458685857
## 8          chain2    -14.0892094958736
## 9          chain3      3.9706692686823
## 10         chain4   -0.952853103676081
## 11         chain5    -4.38379640225017
## 12         chain7    -7.42942718257215
## 13        chain11    -19.2502799115492
## 14        chain12     12.5379067508233
## 15        chain13    -13.9757292841679
## 16        chain14       1.025200487434
## 17        chain15     7.54807320683419
## 18        chain16    -2.36948737590514
## 19        chain17     -1.4170547342752
## 20        chain18    -13.4261110148514
## 21        chain19     -76.260601221126
## 22        chain22     9.53312450790367
## 23        chain23    0.397425969270333
## 24        chain24    -15.6698956189974
## 25        chain26     66.4201490102843
## 26        chain27    -1.24061882743597
## 27        chain28     18.5406147809286
## 28        chain29      7.3896709141192
## 29        chain30    -17.7907512630858
## 30        chain40     8.67948024079653