<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/cross-fitting-double-machine-learning-estimator/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://insightr.wordpress.com/2017/06/28/cross-fitting-double-machine-learning-estimator/"> R – insightR</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<div><strong>By Gabriel Vasconcelos</strong></div>
<h2>Motivation</h2>
<p>In a <a rel="nofollow" href="https://insightr.wordpress.com/2017/05/12/problems-of-causal-inference-after-selection-of-controls/#more-184" target="_blank">late post</a> I talked about inference after model selection showing that a simple double selection procedure is enough to solve the problem. In this post I’m going to talk about a generalization of the double selection for any Machine Learning (ML) method described by Chernozhukov et al. (2016).</p>
<p><span></span></p>
<p>Suppose you are in the following framework:</p>
<p><img alt="displaystyle y_i=d_itheta +g_0(boldsymbol{z}_i)+u_i" title="displaystyle y_i=d_itheta +g_0(boldsymbol{z}_i)+u_i" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+y_i%3Dd_i%5Ctheta+%2Bg_0%28%5Cboldsymbol%7Bz%7D_i%29%2Bu_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"><br>
<img alt="displaystyle d_i=m_0(boldsymbol{z}_i)+v_i" title="displaystyle d_i=m_0(boldsymbol{z}_i)+v_i" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+d_i%3Dm_0%28%5Cboldsymbol%7Bz%7D_i%29%2Bv_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>where <img alt="theta" title="theta" src="https://s0.wp.com/latex.php?latex=%5Ctheta&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> is the parameter of interest, <img alt="boldsymbol{z}_i" title="boldsymbol{z}_i" src="https://s0.wp.com/latex.php?latex=%5Cboldsymbol%7Bz%7D_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> is a set of control variables and <img alt="u_i" title="u_i" src="https://s0.wp.com/latex.php?latex=u_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and <img alt="v_i" title="v_i" src="https://s0.wp.com/latex.php?latex=v_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> are error terms. The functions <img alt="m_0" title="m_0" src="https://s0.wp.com/latex.php?latex=m_0&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and <img alt="g_0" title="g_0" src="https://s0.wp.com/latex.php?latex=g_0&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> are very generic and possibly non-linear.</p>
<p>How can we estimate <img alt="theta" title="theta" src="https://s0.wp.com/latex.php?latex=%5Ctheta&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">? The most naive way (and obviously wrong) is to simple estimate a regression using only <img alt="d_i" title="d_i" src="https://s0.wp.com/latex.php?latex=d_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> to explain <img alt="y_i" title="y_i" src="https://s0.wp.com/latex.php?latex=y_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. Another way is to try something similar to the double selection, which is referred by Chernozhukov et al. (2016) as Double Machine Learning (DML):</p>
<ul>
<li>(1): Estimate <img alt="d_i=hat{m}_0(boldsymbol{z}_i)+hat{v}_i" title="d_i=hat{m}_0(boldsymbol{z}_i)+hat{v}_i" src="https://s0.wp.com/latex.php?latex=d_i%3D%5Chat%7Bm%7D_0%28%5Cboldsymbol%7Bz%7D_i%29%2B%5Chat%7Bv%7D_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">,</li>
<li>(2): Estimate <img alt="y_i=hat{g}_0(boldsymbol{z}_i)+hat{u}_i" title="y_i=hat{g}_0(boldsymbol{z}_i)+hat{u}_i" src="https://s0.wp.com/latex.php?latex=y_i%3D%5Chat%7Bg%7D_0%28%5Cboldsymbol%7Bz%7D_i%29%2B%5Chat%7Bu%7D_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> without using <img alt="d_i" title="d_i" src="https://s0.wp.com/latex.php?latex=d_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">,</li>
<li>(3): Estimate <img alt="hat{theta}=(sum_{i=1}^N hat{v}_id_i)^{-1}sum_{i=1}^N hat{v}_i (y_i-hat{g}_0(boldsymbol{z}_i))" title="hat{theta}=(sum_{i=1}^N hat{v}_id_i)^{-1}sum_{i=1}^N hat{v}_i (y_i-hat{g}_0(boldsymbol{z}_i))" src="https://s0.wp.com/latex.php?latex=%5Chat%7B%5Ctheta%7D%3D%28%5Csum_%7Bi%3D1%7D%5EN+%5Chat%7Bv%7D_id_i%29%5E%7B-1%7D%5Csum_%7Bi%3D1%7D%5EN+%5Chat%7Bv%7D_i+%28y_i-%5Chat%7Bg%7D_0%28%5Cboldsymbol%7Bz%7D_i%29%29&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">.</li>
</ul>
<p>However, in this case the procedure above will leave a term on the asymptotic distribution of <img alt="hat{theta}" title="hat{theta}" src="https://s0.wp.com/latex.php?latex=%5Chat%7B%5Ctheta%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> that will cause bias. This problem may be solved with sample splitting. Suppose we randomly split our <img alt="N" title="N" src="https://s0.wp.com/latex.php?latex=N&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> observations in two samples of size <img alt="n=N/2" title="n=N/2" src="https://s0.wp.com/latex.php?latex=n%3DN%2F2&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. The fist sample will be indexed by <img alt="I" title="I" src="https://s0.wp.com/latex.php?latex=I&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and the auxiliary second sample will be indexed by <img alt="I^c" title="I^c" src="https://s0.wp.com/latex.php?latex=I%5Ec&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. We are going to estimate the first two steps above in the auxiliary sample <img alt="I^c" title="I^c" src="https://s0.wp.com/latex.php?latex=I%5Ec&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and the third step will be done into sample <img alt="I" title="I" src="https://s0.wp.com/latex.php?latex=I&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. Therefore:</p>
<p><img alt="displaystyle hat{theta}=left(sum_{i=in I} hat{v}_id_i right)^{-1}sum_{iin I} hat{v}_i (y_i-hat{g}_0(boldsymbol{z}_i))" title="displaystyle hat{theta}=left(sum_{i=in I} hat{v}_id_i right)^{-1}sum_{iin I} hat{v}_i (y_i-hat{g}_0(boldsymbol{z}_i))" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+%5Chat%7B%5Ctheta%7D%3D%5Cleft%28%5Csum_%7Bi%3D%5Cin+I%7D+%5Chat%7Bv%7D_id_i+%5Cright%29%5E%7B-1%7D%5Csum_%7Bi%5Cin+I%7D+%5Chat%7Bv%7D_i+%28y_i-%5Chat%7Bg%7D_0%28%5Cboldsymbol%7Bz%7D_i%29%29&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>The estimator above is unbiased. However, you are probably wondering about efficiency loss because <img alt="hat{theta}" title="hat{theta}" src="https://s0.wp.com/latex.php?latex=%5Chat%7B%5Ctheta%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> was estimated using only half of the sample. To solve this problem we must now do the opposite: first we estimate steps 1 and 2 in the sample <img alt="I" title="I" src="https://s0.wp.com/latex.php?latex=I&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and then we estimate <img alt="theta" title="theta" src="https://s0.wp.com/latex.php?latex=%5Ctheta&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> in the sample <img alt="I^c" title="I^c" src="https://s0.wp.com/latex.php?latex=I%5Ec&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. As a result, we will have <img alt="hat{theta}(I^c,I)" title="hat{theta}(I^c,I)" src="https://s0.wp.com/latex.php?latex=%5Chat%7B%5Ctheta%7D%28I%5Ec%2CI%29&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and <img alt="hat{theta}(I,I^c)" title="hat{theta}(I,I^c)" src="https://s0.wp.com/latex.php?latex=%5Chat%7B%5Ctheta%7D%28I%2CI%5Ec%29&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">. The cross-fitting estimator will be:</p>
<p><img alt="displaystyle hat{theta}_{CF}=frac{hat{theta}(I^c,I)+hat{theta}(I,I^c)}{2}" title="displaystyle hat{theta}_{CF}=frac{hat{theta}(I^c,I)+hat{theta}(I,I^c)}{2}" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+%5Chat%7B%5Ctheta%7D_%7BCF%7D%3D%5Cfrac%7B%5Chat%7B%5Ctheta%7D%28I%5Ec%2CI%29%2B%5Chat%7B%5Ctheta%7D%28I%2CI%5Ec%29%7D%7B2%7D&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p>which is a simple average between the two <img alt="theta" title="theta" src="https://s0.wp.com/latex.php?latex=%5Ctheta&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">s.</p>
<h2>Example</h2>
<p>The best way to illustrate this topic is using simulation. I am going to generate data from the following model:</p>
<p><img alt="displaystyle y_i=theta d_i + cos^2(boldsymbol{z}_i' b) + u_i" title="displaystyle y_i=theta d_i + cos^2(boldsymbol{z}_i' b) + u_i" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+y_i%3D%5Ctheta+d_i+%2B+cos%5E2%28%5Cboldsymbol%7Bz%7D_i%27+b%29+%2B+u_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<p><img alt="displaystyle d_i = sin(boldsymbol{z}_i'b)+cos(boldsymbol{z}_i'b)+v_i" title="displaystyle d_i = sin(boldsymbol{z}_i'b)+cos(boldsymbol{z}_i'b)+v_i" src="https://s0.wp.com/latex.php?latex=%5Cdisplaystyle+d_i+%3D+sin%28%5Cboldsymbol%7Bz%7D_i%27b%29%2Bcos%28%5Cboldsymbol%7Bz%7D_i%27b%29%2Bv_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"></p>
<ul>
<li>The number of observations and the number of simulations will be 500,</li>
<li>The number of variables in <img alt="boldsymbol{z}_i" title="boldsymbol{z}_i" src="https://s0.wp.com/latex.php?latex=%5Cboldsymbol%7Bz%7D_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> will be <img alt="K=10" title="K=10" src="https://s0.wp.com/latex.php?latex=K%3D10&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">, generated from a multivariate normal distribution,</li>
<li><img alt="theta=0.5" title="theta=0.5" src="https://s0.wp.com/latex.php?latex=%5Ctheta%3D0.5&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">,</li>
<li><img alt="b_k=frac{1}{k},~~~k=1,dots,K" title="b_k=frac{1}{k},~~~k=1,dots,K" src="https://s0.wp.com/latex.php?latex=b_k%3D%5Cfrac%7B1%7D%7Bk%7D%2C~~~k%3D1%2C%5Cdots%2CK&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">,</li>
<li><img alt="u_i" title="u_i" src="https://s0.wp.com/latex.php?latex=u_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and <img alt="v_i" title="v_i" src="https://s0.wp.com/latex.php?latex=v_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> are normal with mean 0 and variance 1.</li>
</ul>
<pre>library(clusterGeneration)
library(mvtnorm)
library(randomForest)

set.seed(123) # = Seed for Replication = #
N=500 # = Number of observations = #
k=10 # = Number of variables in z = #
theta=0.5
b=1/(1:k)

# = Generate covariance matrix of z = #
sigma=genPositiveDefMat(k,&quot;unifcorrmat&quot;)$Sigma
sigma=cov2cor(sigma)
</pre>
<p>The ML model we are going to use to estimate steps 1 and 2 is the Random Forest. The simulation will estimate the simple OLS using only <img alt="d_i" title="d_i" src="https://s0.wp.com/latex.php?latex=d_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> to explain <img alt="y_i" title="y_i" src="https://s0.wp.com/latex.php?latex=y_i&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex">, the naive DML without sample splitting and the Cross-fitting DML. The 500 simulations may take a few minutes.</p>
<pre>set.seed(123)
M=500 # = Number of Simumations = #

# = Matrix to store results = #
thetahat=matrix(NA,M,3)
colnames(thetahat)=c(&quot;OLS&quot;,&quot;Naive DML&quot;,&quot;Cross-fiting DML&quot;)

for(i in 1:M){
  z=rmvnorm(N,sigma=sigma) # = Generate z = #
  g=as.vector(cos(z%*%b)^2) # = Generate the function g = #
  m=as.vector(sin(z%*%b)+cos(z%*%b)) # = Generate the function m = #
  d=m+rnorm(N) # = Generate d = #
  y=theta*d+g+rnorm(N) # = Generate y = #

  # = OLS estimate = #
  OLS=coef(lm(y~d))[2]
  thetahat[i,1]=OLS

  # = Naive DML = #
  # = Compute ghat = #
  model=randomForest(z,y,maxnodes = 20)
  G=predict(model,z)
  # = Compute mhat = #
  modeld=randomForest(z,d,maxnodes = 20)
  M=predict(modeld,z)
  # = compute vhat as the residuals of the second model = #
  V=d-M
  # = Compute DML theta = #
  theta_nv=mean(V*(y-G))/mean(V*d)
  thetahat[i,2]=theta_nv

  # = Cross-fitting DML = #
  # = Split sample = #
  I=sort(sample(1:N,N/2))
  IC=setdiff(1:N,I)
  # = compute ghat on both sample = #
  model1=randomForest(z[IC,],y[IC],maxnodes = 10)
  model2=randomForest(z[I,],y[I], maxnodes = 10)
  G1=predict(model1,z[I,])
  G2=predict(model2,z[IC,])

  # = Compute mhat and vhat on both samples = #
  modeld1=randomForest(z[IC,],d[IC],maxnodes = 10)
  modeld2=randomForest(z[I,],d[I],maxnodes = 10)
  M1=predict(modeld1,z[I,])
  M2=predict(modeld2,z[IC,])
  V1=d[I]-M1
  V2=d[IC]-M2

  # = Compute Cross-Fitting DML theta
  theta1=mean(V1*(y[I]-G1))/mean(V1*d[I])
  theta2=mean(V2*(y[IC]-G2))/mean(V2*d[IC])
  theta_cf=mean(c(theta1,theta2))
  thetahat[i,3]=theta_cf

}

colMeans(thetahat) # = check the average theta for all models = #
</pre>
<pre>##              OLS        Naive DML Cross-fiting DML
##        0.5465718        0.4155583        0.5065751
</pre>
<pre># = plot distributions = #
plot(density(thetahat[,1]),xlim=c(0.3,0.7),ylim=c(0,14))
lines(density(thetahat[,2]),col=2)
lines(density(thetahat[,3]),col=4)
abline(v=0.5,lty=2,col=3)
legend(&quot;topleft&quot;,legend=c(&quot;OLS&quot;,&quot;Naive DML&quot;,&quot;Cross-fiting DML&quot;),col=c(1,2,4),lty=1,cex=0.7,seg.len = 0.7,bty=&quot;n&quot;)
</pre>
<p><img alt="plot of chunk unnamed-chunk-8" src="https://i0.wp.com/i.imgur.com/lqf7qIT.png"></p>
<p>As you can see, the only unbiased distribution is the Cross-Fitting DML. However, the model used to estimate the functions <img alt="m_0" title="m_0" src="https://s0.wp.com/latex.php?latex=m_0&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> and <img alt="g_0" title="g_0" src="https://s0.wp.com/latex.php?latex=g_0&bg=ffffff&%23038;fg=333333&%23038;s=0" class="latex"> must be able to capture the relevant information for the data. In the linear case you may use the LASSO and achieve a result similar to the double selection. Finally, what we did here was a 2-fold Cross-Fitting, but you can also do a k-fold Cross-Fitting just like you do a k-fold cross-validation. The size of k does not matter asymptotically, but in small samples the results may change.</p>
<h2>References</h2>
<p>Chernozhukov, V., Chetverikov, D., Demirer, M., Duflo, E., &amp; Hansen, C. (2016). Double machine learning for treatment and causal parameters. arXiv preprint arXiv:1608.00060.</p>
<p> <a rel="nofollow" href="http://feeds.wordpress.com/1.0/gocomments/insightr.wordpress.com/429/" target="_blank"><img border="0" alt="" src="http://feeds.wordpress.com/1.0/comments/insightr.wordpress.com/429/"></a> <img border="0" alt="" src="https://i1.wp.com/pixel.wp.com/b.gif?resize=1%2C1&ssl=1" data-recalc-dims="1"></p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://insightr.wordpress.com/2017/06/28/cross-fitting-double-machine-learning-estimator/"> R – insightR</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img width="1" alt="" src="http://feeds.feedburner.com/~r/RBloggers/~4/57UCgkuNQwk" height="1"></a></div>