This paper provides estimation and inference methods for a large number of heterogeneous treatment effects in the presence of an even larger number of controls and unobserved unit heterogeneity. In our main example, the vector of heterogeneous treatments is generated by interacting the base treatment variable with a subset of controls. We first estimate the unit-specific expectation functions of the outcome and each treatment interaction conditional on controls and take the residuals. Second, we report the Lasso (L1-regularized least squares) estimate of the heterogeneous treatment effect parameter, regressing the outcome residual on the vector of treatment ones. We debias the Lasso estimator to conduct simultaneous inference on the target parameter by Gaussian bootstrap. We account for the unobserved unit heterogeneity by projecting it onto the time-invariant covariates, following the correlated random effects approach of Mundlak (1978) and Chamberlain (1982). We demonstrate our method by estimating price elasticities of groceries based on scanner data.