The factor of two is coming from the place you identified.
Think about throwing out that factor of two, so you're considering only the bottom hemisphere. When you make your Gaussian shell and have it enclose charge in the bottom hemisphere only, the charge is no longer uniformly distributed inside your Gaussian shell. Thus, the electric field created by the charge you're considering is not the same at all parts of the shell, so you can't find the magnitude of the electric field in the way you described. That only works when the charge distribution has some sort of symmetry you're exploiting. You'd have to do a difficult integral instead.
However, if you don't throw out that factor of two, you're simply finding the electric field inside the shell. Suppose you carry out the rest of your calculation. Then you've found the net force in the z-direction in the north half of the sphere. However, the north half cannot exert any net force on itself, so this entire net force must be the same as the net force from the southern hemisphere.
So you're including all the charge when you make your Gaussian surface because you need to find the true electric field in the shell. The true electric field, when integrated, gives you the net force, which by basic mechanics arguments must be due to the southern hemisphere.
First of all, let's see what Gauss's divergence theorem tells: the outward flux of a vector field through a closed surface is equal to the volume integral of the divergence over the region inside the surface. Intuitively, it states that the sum of all sources minus the sum of all sinks gives the net flow out of a region.
Now, let's look at the Gauss's law in electrostatics:
In differential form, it reads
$$\nabla\cdot\vec{E}=\frac{\rho_{enc}}{\epsilon_0}$$
This means the net outward flux of the electric field lines normal to the surface enclosing the charge is equal to the net charge enclosed by the surface. On integrating the above equation over a spherical volume enclosing the charge,
$$\int_V \nabla\cdot\vec{E}d\tau'=\frac{1}{\epsilon_0}\int_V \rho_{enc}d\tau'$$
By Gauss's divergence theorem, this volume integral of $\vec{E}$ is equal to the outward flux of $\vec{E}$ throgh a closed surface enclosing the charge:
$$\int_V \nabla\cdot\vec{E}d\tau'=\int_{\sigma}\vec{E}\cdot d\vec{\sigma}$$
Hence
$$\int_{\sigma}\vec{E}\cdot d\vec{\sigma}=\frac{1}{\epsilon_0}\int_V \rho_{enc}d\tau'=\frac{q_{enc}}{\epsilon_0}$$
where we have assumed that the volume charge density is continuous and constant. This is Gauss's law in integral form.
So, to use Gauss's law, you should choose the integrating region to be a surface that encloses the charge.
Now, let's look at your problem.
To find the electric field at some point outside the sphere of radius $R$:
We have
$$\int_{\sigma}\vec{E}\cdot d\vec{\sigma}=\frac{q_{enc}}{\epsilon_0}$$
where the integration is over a Gaussian spherical surface enclosing the charged sphere of radius $r$ such that $r>R$ Since the electric field is symmetrical about a spherical surface, we can take it out of the integral. Also the electric field pointing outside the surface is in the same direction as area vector of the spherical surface. Taking $q_{enc}=q$, the total charge enclosed by the charged sphere:
$$E\int_{\sigma}d\sigma=E.4\pi r^2=\frac{q}{\epsilon_0}$$
$$E=\frac{1}{4\pi\epsilon_0}\frac{q}{r^2}$$
In particular, at the surface ($r=R$),
$$E=\frac{1}{4\pi\epsilon_0}\frac{q}{R^2}$$
To find the electric field at some point inside the sphere of radius $R$:
Here our Gaussian sphere is inside the charged sphere. i.e., $r<R$. Since the Gaussian surface is inside the charged sphere, our Gaussian surface do not enclose the entire charge distribution, but only less. Hence the electric field at any point inside the charged sphere will be less than that at the surface.
First we need to find out what is $q_{enc}$.
$$q_{enc}=\int_0^r \rho d\tau'=q\left(\frac{r}{R}\right)^3$$.
Hence the electric field inside the sphere is
$$E.4\pi r^2=\frac{1}{\epsilon_0}q\left(\frac{r}{R}\right)^3$$
$$E=\frac{1}{4\pi\epsilon_0}\frac{qr}{R^3} $$
So,
$$
\vec{E}=
\begin{cases}
\Large{\frac{1}{4\pi\epsilon_0}\frac{q}{r^2}\hat{r}}, &\text{if $r>R$; outside the charged sphere}\\[2ex]
\Large{\frac{1}{4\pi\epsilon_0}\frac{q}{R^2}\hat{r}}, &\text{if $r=R$; on the charged sphere}\\[2ex]
\Large{\frac{1}{4\pi\epsilon_0}\frac{qr}{R^3}\hat{r}}, &\text{if $r<R$; inside the charged sphere}
\end{cases}
$$
Best Answer
I think you forgot to account for $\mathbf{\hat{z}}$
Let $\mathbf{\hat{r}}$ be the normal to the surface of our sphere. If you take the route of integrating the electric field over the surface of the sphere that contains the charge, then you will be evaluating the following quantity.
$z^2 \mathbf{\hat{z}}.\mathbf{\hat{r}}=z^2 \mathbf{\hat{z}}.\mathbf{r}/r=z^2\, z/r=z^3/r$
So you will be integrating $z^3$ over the surface of the sphere centered on the origin. Since $z^3$ is an odd function the integral will vanish.