Find the Var$(X)$ given that $$m_X(u)=\Big(\frac{p}{1-(1-p)e^u)}\Big)^r \ \ \ \ \ \ \ u<\text{ln}((1-p)^{-1})$$
I have found $\mathbb{E}(X)$ to be $$\mathbb{E}(X)=m'_X(u)$$ $$\mathbb{E}(X)=r\Bigg(\frac{p}{(1-(1-p)e^u)}\Bigg)^{r-1}\times\frac{p(e^u-pe^u)}{(1-(1-p)e^u)^2}$$ substituting $u=0$ $$\mathbb{E}(X)=r\Bigg(\frac{p}{p}\Bigg)^{r-1}\times\frac{p(1-p)}{p^2}=\frac{r(1-p)}{p}$$
My question is, how do I find $\mathbb{E}(X^2)$? Differentiating $\mathbb{E}(X)$ again seems like it's impossible.
Sounds like the issue is with calculating the second derivative, so here are some tips. When I have to do a tedious derivative, I like to fold up the notation a little to make things easier. So in this case, I would define
$$ f(u)\equiv 1-(1-p)e^u, $$
then the MGF is
$$ m_X(u)=\bigg(\frac{p}{f(u)}\bigg)^r. $$
Using this notation, we have that
$$ m'_X(u)=-\frac{rp^rf'(u)}{\big[f(u)\big]^{r+1}}. $$
This is the tedious part. Using the quotient rule, the second derivative comes out to
$$ m''_X(u)=\frac{rp^r\big[(r+1)[f'(u)]^2-f(u)f''(u)\big]}{\big[f(u)\big]^{r+2}} $$
You can then easily evaluate $f(0)$, $f'(0)$, etc. to get the desired result.
Another way to do this is to observe that there is a relationship between the moment generating function (MGF) and the probability generating function (PGF); namely, $$m_X(\log u) = \operatorname{E}[e^{X \log u}] = \operatorname{E}[e^{\log u^X}] = \operatorname{E}[u^X] = P_X(u).$$ Furthermore, the PGF has the property that $$\operatorname{E}[X(X-1)\ldots(X-n+1)] = \frac{d^n}{du^n}\left[P_X(u)\right]_{u=1}.$$ So for the negative binomial case, $$P_X(u) = \left(\frac{p}{1-(1-p)u}\right)^r,$$ and the first derivative is $$\frac{dP}{du} = p^r (-r)(-(1-p))(1-(1-p)u)^{-r-1} = \frac{r(1-p)}{p} \left(\frac{p}{1-(1-p)u}\right)^{r+1}.$$ Evaluating at $u = 1$ gives $$\operatorname{E}[X] = \frac{r(1-p)}{p}$$ as desired. But now note $$\frac{dP}{du} = \operatorname{E}[X] P_{X^*}(u),$$ where $X^* \sim \operatorname{NegBinomial}(r+1,p)$. So the second derivative is trivial: $$\frac{d^2 P}{du^2} = \operatorname{E}[X] \operatorname{E}[X^*] P_{X^{**}}(u),$$ where $X^{**} \sim \operatorname{NegBinomial}(r+2,p)$, and the pattern continues, letting us conclude in general that $$\frac{d^n}{du^n}\left[P_X(u)\right]_{u=1} = \prod_{k=0}^{n-1} (r+k) \frac{1-p}{p} = \left(\frac{1-p}{p}\right)^n \frac{(r+n-1)!}{(r-1)!},$$ which is rather nice since now $$\operatorname{E}\left[\binom{X}{n}\right] = \binom{r+n-1}{n} \left(\frac{1-p}{p}\right)^n.$$ This is quite above and beyond the original request, but we can now easily compute the variance: $$\begin{align*} \operatorname{Var}[X] &= \operatorname{E}[X^2]-\operatorname{E}[X]^2 \\ &= \operatorname{E}[X(X-1)+X] - \operatorname{E}[X]^2 \\ &= \operatorname{E}[X(X-1)]+\operatorname{E}[X] - \operatorname{E}[X]^2 \\ &= \operatorname{E}[X](\operatorname{E}[X^*] + 1 - \operatorname{E}[X]) \\ &= \frac{r(1-p)}{p} \left( \frac{(r+1)(1-p)}{p} + 1 - \frac{r(1-p)}{p}\right) \\ &= \frac{r(1-p)}{p^2} ((r+1)(1-p) - r(1-p) + p) \\ &= \frac{r(1-p)}{p^2} (1-p + p) \\ &= \frac{r(1-p)}{p^2}. \end{align*}$$
Let $q:=1-p$
$\begin{align}\mathsf M'_X(u) & = \mathcal D_u\left(\dfrac{p}{1-(1-p)e^u}\right)^r\\[1ex] & =p^r \mathcal D_u(1-qe^u)^{-r}\\[1ex] &= -rp^r (1-qe^u)^{-r-1}\mathcal D_u (1-qe^u)\\[1ex] &= rqp^r (1-qe^u)^{-r-1}e^u\\[2ex]\mathsf M''_X(u) &= \mathcal D_u\left(rqp^r (1-qe^u)^{-r-1}e^u\right)\\ & =rqp^r\left(e^u\mathcal D_u(1-qe^u)^{-r-1} + (1-qe^u)^{-r-1}\mathcal D_u e^u\right)\\ & = rqp^r\left((r+1)e^u(1-qe^u)^{-r-2}qe^u + (1-qe^u)^{-r-1}e^u\right)\\&=rqp^r(1-qe^u)^{-r-2}\left((r+1)qe^{2u} + e^u-qe^{2u}\right)\\&=rqp^r(1-qe^u)^{-r-2}\left(rqe^{2u}+e^u\right)\\[2ex] \mathsf M'_X(0) & = rqp^{-1}\\&=r(1-p)p^{-1}\\[1ex] \mathsf M''_{X}(0) & = (r^2q^2+rq)\,p^{-2} \\[1ex] \mathsf M''_X(0)-(\mathsf M'_X(x))^2 &=rqp^2 \\&=r(1-p)p^{-2} \end{align}$
I wanted to state that finding the variance of $X$, given the MGF of $X$, is usually much easier done with the cumulant generating function (CGF).
Suppose $M_X$ is the MGF of $X$. Then the CGF is given by $\varphi_X = \ln M_X$.
One of the very nice things about $\varphi_X$ is that $\left. \varphi^{\prime\prime}_X(t) \right|_{t = 0} = \text{Var}(X)$, thus eliminating the need to find $\mathbb{E}[X]$ and $\mathbb{E}[X^2]$ separately.
Given $$M_X(u) = \left[ \dfrac{p}{1-(1-p)e^u}\right]^r \implies\varphi_X(u)=\ln M_X(u)=r\ln\left[ \dfrac{p}{1-(1-p)e^u}\right]$$ where we have used the property $\ln(a^b) = b\ln(a)$.
The first derivative of $\varphi_X$ is given by $$\varphi_X^{\prime}(u) = r\dfrac{1}{p/[1-(1-p)e^u]} \cdot p \cdot [1-(1-p)e^u]^{-2}[(1-p)e^{u}] =\dfrac{r(1-p)e^u}{1-(1-p)e^u}$$ [It's worth noting that $\left.\varphi^{\prime}_X(u)\right|_{u = 0} = \mathbb{E}[X]$.]
To make finding the second derivative easier, divide both the numerator and denominator by $e^u$, obtaining $$\varphi^{\prime}_X(u)=\dfrac{r(1-p)}{e^{-u}-(1-p)}$$ The second derivative is then $$\varphi^{\prime\prime}_X(u)=r(1-p)(-1)[e^{-u}-(1-p)]^{-2}(-e^{-u})=\dfrac{re^{-u}(1-p)}{[e^{-u}-(1-p)]^2}$$ and evaluating this at $u = 0$ gives $$\left.\varphi^{\prime\prime}_X(u)\right|_{u=0}=\dfrac{re^{-0}(1-p)}{[e^{-0}-(1-p)]^2} = \dfrac{r(1-p)}{p^2} = \text{Var}(X)$$ as desired.
