Operator $\rho$ is not a tensor product, it's a sum of tensor products
$$
p_1|1\rangle\langle 1| \otimes \rho_1 + p_2|2\rangle\langle 2| \otimes \rho_2 + \dots + p_d|d\rangle\langle d| \otimes \rho_d.
$$
This is not the same as
$$
\big(\sum_ip_i|i\rangle\langle i|\big) \otimes \big(\sum_i\rho_i\big),
$$
so your expansion isn't correct.
Also in general $S(A+B)\neq S(A)+S(B)$, but in this situation the supports of $|i \rangle\langle i|\otimes \rho_i$ and $|j \rangle\langle j|\otimes \rho_j$ are orthogonal, so we can write
$$
S(\rho) = S(p_1|1\rangle\langle 1| \otimes \rho_1) + \dots + S(p_d|d\rangle\langle d| \otimes \rho_d)
$$
Here $p_i|1\rangle\langle 1| \otimes \rho_i$ is not a density matrix because it's scaled, i.e. its trace equals $p_i<1$, so technically $S$ is not defined. But for such matrices we also can define expression $S(M) = -\sum_i \lambda_i\text{ln}\lambda_i$, where $\lambda_i$ are eigenvalues of $M$. It's easy to check that for $c>0$ and density matrix $\rho$ we have $S(c\rho) = cS(\rho) - c\text{ln}c$. So we can write
$$
S(p_i|i\rangle\langle i| \otimes \rho_i) = p_iS(|i\rangle\langle i| \otimes \rho_i)-p_i\text{ln}p_i =
$$
$$
= p_i\big(S(|i\rangle\langle i|) + S(\rho_i)\big)-p_i\text{ln}p_i = p_iS(\rho_i)-p_i\text{ln}p_i
$$
After the summation we will have
$$
S(\rho) = \sum _i p_iS(\rho_i) - \sum_ip_i\text{ln}p_i = \sum _i p_iS(\rho_i) + S(\overline{p})
$$