OLS regression and dummy variables: fitted values of the subsample equal actual values? Given $X=(d\; X_1)$, we want to prove that
$$(d'd)^{-1}d'P_{[d X_1]}y=(d'd)^{-1}d'y,$$ where $P_{[X]}y=\hat{y}.$ That is, restricting the regression to the subsample for which $d_i=1$, we have that the OLS fitted values equal the actual values.
I have a problem proving this.
(Attempted) Proof:
Given that $P_{[X]}=P_{[d]}+P_{[M_{[d]}X_1]}$, we have:
$$(d'd)^{-1}d'P_{[d X_1]}y=(d'd)^{-1}d'(P_{[d]}+P_{[M_{[d]}X_1]})y
\\=(d'd)^{-1}d'P_{[d]}y + (d'd)^{-1}d'P_{[M_{[d]}X_1]}y
\\= (d'd)^{-1}d'd(d'd)^{-1}d'y + (d'd)^{-1}d'P_{[M_{[d]}X_1]}y
\\= (d'd)^{-1}d'y + (d'd)^{-1}d'P_{[M_{[d]}X_1]}y$$
and here is where I am at a loss. It appears that $(d'd)^{-1}d'P_{[M_{[d]}X_1]}y$ should equal 0, in this way we would get the result at the beginning. The only way for this to happen, apart from the trivial case where $y=0$ or $d=0$, is that $P_{[M_{[d]}X_1]}y=0$. That is, only if $M_{[d]}X_1$ and $y$ where orthogonal.
Is this line of reasoning correct? And, in any case, what would be the meaning of this supposed orthogonality?

Since $d_i\in\{0,1\}$, $d'd=n_1$, where $n_1$ denotes the number of observations for which $d_i=1$. Also,
$$
d'y=\sum_{d_i=1}y_i
$$
so that
$$
(d'd)^{-1}d'y=\frac{1}{n_1}\sum_{d_i=1}y_i
$$
is the mean of the fitted values of the regression. Along the same lines, the lhs of your first display is the mean of the fitted values of the full regression for which $d_i=1$.
Hence, you are to prove that two numbers are identical, not two vectors (which, as Dave wrote, would be too good to be true...).
It is therefore also enough to prove
$$
d'P_{[d X_1]}y=d'y.
$$
Here, $d'P_{[d X_1]}=(P_{[d X_1]}'d)'=(P_{[d X_1]}d)'$ by symmetry of projection matrices.
These are the (transposed) fitted values of a regression of $d$ on $d$ and $X_1$. Since this regression explains $d$ by, inter alia, itself, we have (see Hayashi Econometrics Seemingly Unrelated Regressions (SUR) Eq 4.5.13'-15' for a formal argument)
$P_{[d X_1]}d=d$ and hence
$$d'P_{[d X_1]}y=d'y.$$
Your proof also seems correct as far as it goes. At the end, you have
$$
d'P_{[M_{[d]}X_1]}
$$
Here, $P_{[M_{[d]}X_1]}$ starts with $M_{[d]}$, so that $d'M_{[d]}=0$, and we do not need that $P_{[M_{[d]}X_1]}y=0$.
Numerical verification:
n <- 20
d <- sample(c(0,1), n, replace=T)
x1 <- rnorm(n)
y <- rnorm(n)

fullreg <- lm(y~cbind(d, x1)-1)

mean(y[d==1])
mean(fitted(fullreg)[d==1])

# these are also the coefficient estimates from the following regressions:
(dreg <- summary(lm(y~d-1)))
(dreg2 <- summary(lm(fitted(fullreg)~d-1)))

