Introduction

This is my attempte to replication the table III of Maxim Engers, Monica Hartmann, and Steven Stern (2009) by using Fortran. This paper investigates what factors influence household miles decisions and whether these choices are consistent with the observed pricing paths of these assets. The main contributions of this paper are:

  • The results in this paper indicates that the structural estimates of the household mileage decision better explain the observed pattern of price decline over a vehicle’s life than the OLS model.
  • The structural model they constructed allows decomosing the age effect of annual mileage on prices into three parts: a direct aging effect, a household–car portfolio effect, and a household demographics effect.
  • The superiority of the structural model suggests that one cannot estimate the effect age has on annual miles decisions independent of household characteristics and household–car portfolio choices.
  • Their model of vehicle utilization and used car prices differs from pervious studies.

The following sections provide a very brief overview of the model in this paper, the data source and the corresponding Fortran codes for how it is estimated.

1. Model Description

This paper presents a framework which enables one to determine whether variantion in vehicle usage helps explian the observed price declines for automobiles. In this framework, they first model the household vehicle annual miles decision and then establish the empirical relationship between annual miles and used car prices. To estimate the household-vehicle annual miles decision, they estimate two models: a basic OLS model and a structural model. The OLS model assumes a linear relationship between household characteristics and annual miles. In contrast, the structural model allows for nonlinear influence of changes in car portfolios and household demographic characteristics on how a vehicle’s age affects annual miles.

1.1 Dynamic model of annual miles and price

Before introducing the model, first note that throughout the paper, they are interested in the flow of miles driven in one year rather than in the stock of miles previously driven. Suppose that the benefit flow net of costs derived from owning a car of brand \(b\) and age \(t\) is denoted by \(f_{bt}\) which monotonically decreases with \(t\). Suppose that the rate of decline in \(f_{bt}\) is uniform over time (but can vary from brand to brand), so that

\[f_{bt} = \zeta^{f}_b + \xi^f_{bt}\] where \(\zeta^{f}_b > 0\) and \(\xi^f_b t<0\). If \(\beta\) is the per-period discount factor, the value \(W_{bt}\) to any consumer of a car of brand \(b\) of age \(t\) satisfies a simple recursion condition:

\[W_{b t}=\max \left\{f_{b t}+\beta W_{b t+1}, p_{b} t\right\},\]

that is, the car can either be sold at price \(p_{bt}\) or kept for the period, generating current net surplus of \(f_{bt}\). Because \(f_{b0}\) > 0 and f_{bt} decreases at a constant rate, there is a (unique) final period \(t^*_b\) in which the vehicle provides a positive service flow:

\[\exists t_{b}^{*}: p_{b t_{b}^{*}}=W_{b t_{b}^{*}}>0, p_{b t_{b}^{*}+k}=W_{b t_{b}^{*}+k}=0 \forall k>0\]

Let \(\Delta W_{bs}\) denote the change in Wbs when s is increased by one. In all periods \(t \leq t^*_b\), \(W_{b t} = f_{b t}+\beta W_{b t+1}\), and so

\[\Delta W_{b t}=\xi_{b}^{f}+\beta \Delta W_{b t+1} \forall t<t_{b}^{*}\]

By induction on \(t^*_b - t\), it follows that \(\Delta W_{b t} < 0, \forall t<t_{b}^{*}\) so prices decline monotonically until the car becomes worthless. Differentiating gives

\[\frac{\partial \Delta W_{b t+1}}{\partial \xi_{b}^{f}}=1+\beta \frac{\partial \Delta W_{b t+1}}{\partial \xi_{b}^{f}} \forall t<t_{b}^{*}\] and

\[\frac{\partial \Delta W_{b t+1}}{\partial \xi_{b}^{f}}>0 \forall t<t_{b}^{*}\] This equation says that if we compare two different brands then the one whose annual miles decline at a faster rate will also be the one whose price declines faster. The goal of this paper is to see whether this relationship between the annual miles changes and price changes is found in the data we examine.

1.2 Basic OLS model

The specification of the estimated model is

\[\begin{aligned} \log m_{i j}=& \lambda_{00}+\lambda_{10} a_{i j}+\sum_{k=1}^{K} \lambda_{0 k} b_{i j k}+\sum_{k=1}^{K} \lambda_{1 k} b_{i j k} a_{i j} \\ &+\sum_{k=1}^{K} \lambda_{2 k} b_{i j k} \min \left(a_{i j}, 5\right)+\sum_{h=1}^{H} \lambda_{3 h} D_{i h}+e_{i}+\varepsilon_{i j} \\ e_{i} \sim & \text { i.i.d. }\left(0, \sigma_{h}^{2}\right), \varepsilon_{i j} \sim \text { i.i.d. }\left(0, \sigma_{v}^{2}\right) \end{aligned}\]

where \(m_{ij}\) is annual miles of vehicle \(j\) in household \(i\), \(b_{ijk}\) is a dummy for whether vehicle \(j\) is brand \(k\), \(a_{ij}\) is the age of vehicle \(j\), and \(D_i\) is a vector of \(H\) household characteristics.

1.3 Structural Estimation Model

This section presents a structural model that formally incorporates the nonlinear relationship between annual miles and number of drivers. In the model, the annual miles choices are determined by a two-stage process that first determines the number of trips and then determines the choice of car for each trip. The demand for trips is a function of the household’s demographic characteristics and the quality of vehicles it owns. The choice of car to use for each trip is a random utility choice where the value of each mode for the particular trip has a deterministic component plus an i.i.d. Extreme Value error. Here, they models the vehicle choice for each trip and then aggregates over trips to find annual miles shares for each household vehicle. Let \(m_{ij}\) be annual miles of car \(j\) in household \(i\), \(j = 1, ..., J_i\), and \(m_{i\cdot} = \sum^{J_{i}}_{j=1} m_{ij\cdot}\). Let household i have demographic characteristics \(D_i\) and car \(j\) of household \(i\) have characteristics \(C_{ij}\). They can model

\[\begin{array}{l} m_{i j}=g\left(j, C_{i}, u_{i}, \theta\right) m_{i} \\ m_{i}=h\left(C_{i}, D_{i}, u_{i}, \theta\right) \end{array}\] where \(C_i = (C_{i1}, C_{2i}, ..., C_{iJ_i})\), \(u_{i} = (u_{i1},...,u_{iJ_i})\) is a vector of errors, \(u_i \sim i.i.d. N(0, \sigma^2_u I)\), and \(\theta\) is a vector of parameters specified below. The value of car \(j\) can be \[V_{ij} = C_{ij}\alpha + u_{ij}\] and let \(R_i(j)\) be the rank of car \(j\) in household \(i\) according to \(V_i\). Let \(S(\cdot)\) be the inverse of \(R_i(\cdot)\). They model \(h\) semiparametrically using polynomials with imposed monotonicity. In particular, define \(p_0(v) = 1\) and

\[p_{j}(v)=\sum_{m=1}^{M} \delta_{j m} v^{m}\] as an Mth-degree polynomial in \(v\) for \(j > 0\). Let

\[\Delta V_{i S(j)}=\left\{\begin{array}{cl} 0 & \text { if } j=0 \\ V_{i S(j)} & \text { if } j=1 \\ V_{i S(j)}-V_{i S(j-1)} & \text { if } j>1 \end{array}\right.\] where \(\Delta V_{iS(1)}\) is the value of \(i\)’s best choice, and \(\Delta V_{iS(j)}\) is the negative of the marginal reduction in value between the \(j\)th best choice and the \((j - 1)\)th best choice for \(j > 1\). Define \(J^*_i = \min (J_i, 3)\) and

\[\Phi\left[J_{i}^{*}, V_{i}\right]=\sum_{j=0}^{J_{i}^{*}} \sum_{k=j}^{J_{i}^{*}} \kappa_{j k} p_{j}\left(\Delta V_{i S(j)}\right) p_{k}\left(\Delta V_{i S(k)}\right)\]

as a product of polynomials in different combinations of \(\Delta V_{iS(j)}\) and \(\Delta V_{iS(k)}\) with \(\kappa_{00} = 0\) and \(\kappa_{01} = 1\). \(\Phi\left[J_{i}^{*}, V_{i}\right]\) can be thought as a as a semiparametric approximation of a general function in \((\Delta V_{iS(j)})^3_{j=0}\) in the sense that, as the sample size increases, they can increase M. For any function \(f: W \rightarrow R\) where \(W\) is a nonempty subset of \(R^J\) for some \(J\), define

\[\Psi[f](x)=\sup _{y<x}[f(y)]\] \(\Psi\) transforms the function \(f\) into a new function \(\Psi[f]\) that nondecreasing in all of its arguments. They use the \(\Psi[\cdot]\) transformation because annual miles by household should be increasing in the vector of \(V\) values. They model \[\log h\left(C_{i}, D_{i}, u_{i}, \theta\right)=\Psi[\Phi]+D_{i} \gamma \qquad (1) \] as a partial linear equation, flexible in car values and linear in household characteristics. They approximate the \(\Psi\) function by discretizing the \(V\)-space and then checking for monotonicity in \(h\) point by point.

For the function \(g\), they assume it depends on \(C_i\) and \(u_i\) only through the effect of \(C_i\) and \(u_i\) on \(V_i = (V_{i1}, ..., V_i J_i)\), and it depends on \(j\) only through \(V_{ij}\) and \(R_{ij} = R_i(j).\) Let

\[g\left(j, C_{i}, u_{i}, \theta\right)=\frac{\exp \left\{\psi\left(V_{i j}, R_{i j}, \theta\right)\right\}}{\sum_{i=1} \exp \left\{\psi\left(V_{i j^{\prime}}, R_{i j^{\prime}}, \theta\right)\right\}} \qquad (2) \] where

\[\psi\left(V_{i j}, R_{i j}, \theta\right)=\left\{\begin{array}{cl} \sum_{k=R_{i j}}^{J_{i}-1} \omega_{k}\left(V_{i S(k)}-V_{i S(k+1)}\right) & \text { if } R_{i j} \leq J_{i}-1 \\ 0 & \text { if } R_{i j}=J_{i} \end{array}\right. \qquad (3) \]

is a spline function in \(V_{ij}\) with nodes at the values of the ranked \(V_i\)s. This is a weighted logit function that allows for flexibility with respect to how ranking interacts with car value. Note that it is increasing in each element of \(V_i\) and that, if \(w_k = 1 \forall k\), then \(g(j, C_i, u_i, \theta)\) is logit function.

Finally, the functional forms of equations \((1) - (3)\) satisfy the conditions described in Berry et al. (1995) for a unique solution in \(u_i\) to the set of equations

\[\frac{m_{i j}}{m_{i}}=g\left(j, C_{i}, u_{i}, \theta\right), j=1,2, \ldots, J_{i} \qquad (3)\] \[\log m_{i .}=\log h\left(C_{i}, D_{i}, u_{i}, \theta\right) \qquad (4)\]

Define \(\hat{u}_i\) as the value of \(u_i\) that satisfies equations (3) and (4). If we decompose \(\hat{u}_i\) into \(\hat{u}_{i1}\) and \(\hat{u}_{i/1} = (\hat{u}_{i2}, ..., \hat{u}_{i}J_i)\), then, conditional on any value of \(u_{i1} = u_{i1}^0\), there is a unique solution \(u_{i/1}^0\) to equation (3). Furthermore,

\[\frac{\partial u_{i j}^{0}}{\partial u_{i 1}^{0}}=1 \text { for all } j\] because of the usual location identification problem in random utility models. The value of \(\hat{u}_{i1}\) is the value that solves equation (4).

The set of parameters to estimate is \(\theta = (\alpha, \gamma, \delta, \kappa, \omega, \sigma^2_{u})\), where \(\alpha\) is the effect of car characteristics \(C_{ij}\) on car values \(V_{ij}\), \(\gamma\) is the effect of demographic characteristics \(D_i\) on total household annual miles, \(\delta\) is the set of polynomial coefficients, \(\kappa\) is the set of weights on each of the polynomial products used to determine the total value of the household’s car portfolio, \(\omega\) is the set of weights determining how car shares depend upon relative values, and \(\sigma^2_u\) is the variance of the unobserved household heterogeneity. The log-likelihood contribution of family \(i\) is

\[L_{i}=\sum_{j=1}^{J_{i}} \log \frac{1}{\sigma_{u}} \phi\left[\frac{\widehat{u}_{i j}\left(C_{i}, \theta\right)}{\sigma_{u}}\right]-\log |\Gamma(\theta)|\] where \(\phi[\cdot]\) is the standard normal density function and \(\Gamma(\theta)\) is the Jacobian.

2. Data Source Description

There are a number of data sources used in this paper. These are described below:

  • They use an extract of the 1995 National Personal Transportation Survey (1995). The annual miles data for this study come from this survey. The NPTS is a sample of \(47,293\) vehicles owned by \(24,814\) households. This survey includes some characteristics of each household in the sample such as family income, location characteristics (e.g., urban), and the number of drivers disaggregated by age, gender, and work status. It also includes each vehicle owned by the household, its age, brand, and annual miles driven. The Fortran program dataextract.f can be used to replicate the data. The original data and full documentation is available publicly through, for example, ICPSR. A list of variables pulled from the NPTS is described below:
## NPTS Household record:
## Name of Variable in FORTRAN program | Description | FORTRAN format

ivar(1) availability of public transportation 22x,i2
ivar(2) household id  14x,i8
ivar(3) household income 8x,i2
ivar(4) year of interview 44x,i4
ivar(5) availability of public transportations 6x,i2
ivar(6) member age i4
ivar(7) whether member is a driver i2
ivar(8) whether member is male 2x,i2
ivar(9) whether member works 2x,i2
ivar(10) member age i4
ivar(11) whether member is a driver i2
ivar(12) whether member is male 2x,i2
ivar(13) whether member works 2x,i2
ivar(14) member age i4
ivar(15) whether member is a driver i2
ivar(16) whether member is male 2x,i2
ivar(17) whether member works 2x,i2
ivar(18) member age i4
ivar(19) whether member is a driver i2
ivar(20) whether member is male 2x,i2
ivar(21) whether member works 2x,i2
ivar(22) member age i4
ivar(23) whether member is a driver i2
ivar(24) whether member is male 4x,i2
ivar(25) whether member works 2x,i2
ivar(26) member age i4
ivar(27) whether member is a driver i2
ivar(28) whether member is male 2x,i2
ivar(29) whether member works 2x,i2
ivar(30) member age i4
ivar(31) whether member is a driver i2
ivar(32) whether member is male 2x,i2
ivar(33) whether member works 2x,i2
ivar(34) member age i4
ivar(35) whether member is a driver i2
ivar(36) whether member is male 2x,i2
ivar(37) whether member works 2x,i2
ivar(38) member age i4
ivar(39) whether member is a driver i2
ivar(40) whether member is male 2x,i2
ivar(41) whether member works 2x,i2
ivar(42) member age i4
ivar(43) whether member is a driver i2
ivar(44) whether member is male 2x,i2
ivar(45) whether member works 2x,i2
ivar(46) (not used in analysis) 125x,i4 weight household weight f10.2
ivar(47) urban 4x,i2 aivar alphanumeric urban 140x,a1
  • The price data collected from Kelley Blue Book’s Official Guide for Older Cars (May/June - Central/Eastern Edition, 1995-2000) and the 1995 NPTS to run regressions of price on miles driven as a function of age. The price data are generated from sales reports to the National Automobile Dealer’s Association (NADA) provided by member dealers. For each year in the sample period, the data set includes the average price of each brand of car for each relevant manufacturing year. The Kelley Blue Book variables used in our study are described in the ASCII file variables_mileage_drives_prices.txt, and the data are saved in the space delimitated ASCII file data_mileage_drives_prices.prn. The data from data_mileage_drives_prices.prn is manipulated into a form directly used in the OLS regression and can be read from a rectangular ASCII file called olsdata.prn with 26 records. The code necessary to read the data is in the subroutine decompos1 in the Fortran program nptsestimate.f. There are two variables per record; the first is log(avg miles by brand), and the second is log(average price by brand). The order of the brands is in the Fortran statement:
  • The scrapping data used in this study are from R. L. Polk & Co. The scrapping rates are observed reductions in the stock of cars of each brand as a proportion of the stock the year before. This data set includes scrapping rates only for eight car categories and only for a limited number of relevant vintages, so they use nearby substitutes for other brands and use scrapping rates for cars manufactured in 1990. Also, they ignore the first 2 years of the car because observed scrapping rates are actually negative as car dealers sell their remaining stock.

3. Fortran Codes for Extracting Data and Estimation Steps.

Extract data

First we need to extract data from the national personal transportation survey (npts) for 1995.

 parameter(nvars=47,nchars=4,nvehs=9,nhhmoms=9,npmoms=4,nvmoms=6,
     & nmisss=12,nmodels=26,nages=20)
      implicit real*8(a-h,o-z)
      character*1 aivar,aurban(5)
      character*8 ahname(nhhmoms),amiss(nmisss),amodel(nmodels),
     & aodo(13),aorder(10),apname(npmoms),avname(nvmoms)
      dimension agemodkmom(nmodels,nages,3),agemodmom(nmodels,nages,3),
     & agew(3),amfsmh(10,3),dmfsmh(10,20),hhmom(2,nhhmoms,5),
     & ihchar(10,nchars),ihousetr(18),imiss(2,nmisss),ivar(nvars),
     & jveh(nvehs),jvehtr(nvehs),mfszkp(10),mfszvh(10,20),
     & modyrct(nages,200),modyrfrq(55),mveh(30,nvehs),odrfrq(13),
     & odyrmn(nages,2),pmom(2,npmoms,5),rmfsmh(10,10,2),rmiles(10),
     & vmodel(2,nmodels),vmom(2,nvmoms,5)
      data agew/.4d0,.2d0,.1d0/
      data ahname/'weight','hhinc','busavl','pbtravl','urban1','urban2',
     & 'nhchar','ihco','ihcar'/
      data amiss/'hhinc','busavl','pbtravl','urban1','urban2',
     & 'annmiles','mod year','od mon','od year','od read','veh age',
     & 'truck'/
      data amodel/'isuzu','chrysler','dodge','plymouth','ford',
     & 'mercury','buick','chevrolt','oldsmobl','pontiac','saturn',
     & 'luxamer','luxjapan','luxeurop','honda','mitsubi','mazda',
     & 'nissan','subaru','toyota','volkswg','volvo','geo','hyundai',
     & 'other','truck'/
      data aodo/'< 10','10-20','20-30','30-40','40-50','50-60','60-70',
     & '70-80','80-90','90-100','100-150','150-200','>200'/
      data aorder/'first','second','third','fourth','fifth','sixth',
     & 'seventh','eigth','ninth','tenth'/
      data apname/'age','driver','sex','worker'/
      data avname/'annmiles','model','mod year','od month','od year',
     & 'od read'/
      data aurban/'C','R','S','T','U'/
      data ihousetr/25,75,125,175,225,275,325,375,425,475,525,575,625,
     & 675,725,775,900,1200/
  • Initialize
  • Read household data record.
  • Read vehicle record
   9  read(9,121,end=10)(jveh(i),i=1,2),annm,(jveh(i),i=3,nveh)
 121  format(i8,30x,i6,f8.2,6x,i2,2x,i3,14x,i2,i4,82x,2i4,12x,i6)
 
 ## transform data
 ## household id:
      jhouseid=jveh(1)
      ivco=ivco+1
      ifvmiss=0
      
 ##  annualized miles:
      jvehtr(1)=(((jveh(2)*10000)+annm)/100)+1
      if((jvehtr(1).lt.1).or.(jvehtr(1).gt.1151))then
         imiss(1,6)=imiss(1,6)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,6)=imiss(2,6)+1
            endif
         endif
         
  ## make/model:
         jveh3=jveh(3)
      jveh4=jveh(4)
      jveh5=jveh(5)
      call maketr(jvehtrt,jveh3,jveh4,jveh5)
      if(jvehtrt.eq.26)then
         imiss(1,12)=imiss(1,12)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,12)=imiss(2,12)+1
            endif
         endif
      jvehtr(2)=jvehtrt
      
  ## model year:    
      jvehtr(3)=jveh(6)
      if(jvehtr(3).gt.2010)then
         imiss(1,7)=imiss(1,7)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,7)=imiss(2,7)+1
            endif
         endif
      if(jvehtr(3).le.2010)then
         jvehtrt=jvehtr(3)-1954
         if(jvehtr(2).ne.26)modyrfrq(jvehtrt)=modyrfrq(jvehtrt)+1
         iage=ivar(4)+1901-jvehtr(3)
         if(iage.le.0)then
            if(iage.eq.0)iage=1
            if(iage.lt.0)then
               imiss(1,11)=imiss(1,11)+1
               if(ifvmiss.eq.0)then
                  ifvmiss=1
                  imiss(2,11)=imiss(2,11)+1
                  endif
               endif
            endif
         if(iage.gt.20)iage=20
          if(iage.gt.0)then
            iannt=jvehtr(1)
            if(iannt.le.0)iannt=1
            if(iannt.gt.200)iannt=200
            if(jvehtr(2).ne.26)modyrct(iage,iannt)=modyrct(iage,iannt)+1
            imod=jvehtr(2)
            agemodmom(imod,iage,1)=agemodmom(imod,iage,1)+1.d0
            if(jvehtr(1).le.200)ajvehtr=log(dble(float(jvehtr(1)*100)))
            if(jvehtr(1).gt.200)ajvehtr=log(20000.d0)
            agemodmom(imod,iage,2)=agemodmom(imod,iage,2)+ajvehtr
            agemodmom(imod,iage,3)=agemodmom(imod,iage,3)+(ajvehtr*
     &       ajvehtr)
            ib=iage-2
            if(ib.lt.1)ib=1
            ie=iage+2
            if(ie.gt.nages)ie=nages
            do 55 i=ib,ie
               ia=i-iage
               if(ia.lt.0)ia=-ia
               ia=ia+1
               agemodkmom(imod,i,1)=agemodkmom(imod,i,1)+agew(ia)
               agemodkmom(imod,i,2)=agemodkmom(imod,i,2)+(ajvehtr*
     &          agew(ia))
  55           agemodkmom(imod,i,3)=agemodkmom(imod,i,3)+(ajvehtr*
     &          ajvehtr*agew(ia))
            endif
         endif
         
    ##  odometer stuffL
    ##  (4) month
    ##  (5) year 
    ##  (6) reading
         
         jvehtr(4)=jveh(7)
      if(jvehtr(4).gt.12)then
         imiss(1,8)=imiss(1,8)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,8)=imiss(2,8)+1
            endif
         endif
      jvehtr(5)=jveh(8)
      if(jvehtr(5).gt.96)then
         imiss(1,9)=imiss(1,9)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,9)=imiss(2,9)+1
            endif
         endif
      jvehtr(6)=jveh(9)/100
      if(jvehtr(6).gt.9975)then
         imiss(1,10)=imiss(1,10)+1
         if(ifvmiss.eq.0)then
            ifvmiss=1
            imiss(2,10)=imiss(2,10)+1
            endif
         endif
      if(jvehtr(6).le.9975)then
         jvehtrt=jvehtr(6)/100
         do 41 i=1,10
            if(jvehtrt.le.i)then
               odrfrq(i)=odrfrq(i)+1.d0
               goto 40
               endif
  41        continue
         if(jvehtrt.le.15)then
            odrfrq(11)=odrfrq(11)+1.d0
            goto 40
            endif
         if(jvehtrt.le.20)then
            odrfrq(12)=odrfrq(12)+1.d0
            goto 40
            endif
         odrfrq(13)=odrfrq(13)+1.d0
  40     continue
         if(iage.gt.0)then
            odyrmn(iage,1)=odyrmn(iage,1)+1.d0
            odyrmn(iage,2)=odyrmn(iage,2)+dble(float(jvehtr(6)))
            endif
         endif
  • chech for match.

Step 1: compute the distribution of age derivatives disaggregated by brand and age.

  • Compute derivatives for each observation.
  • Compute errors.
  • compute mileage for each car and household characteristics term.
  • compute total miles for household.
  • compute proportions and mileage for each car and derivatives
  • compute mileage for each car.
  • update derivative moments.
  • adjust moments and output.

Step 2: estimates the average age effect on average mileage disaggregated by number of household cars.

  • initialize.
  • check if the j’th household has the correct number of cars.
  • update relevant moment matrices.

Step 3: computes the covariance matrix of the estimates.

  • Initialize.
  • update with each observation
  • compute log likelihood contribution at estimate.
  • deviate each parameter and compute numerical derivative.
  • update mean and outer product and adjust moments.
  • output average log likelihood contribution derivatives.

4.Results

  • Currently I cannot get the results using the codes above.