Solution (source code)

= Solution

Differentiate the normalization of the <exponential family> density with respect to $\theta_i$, assuming the regularity that permits differentiation under its <integral>. The density <derivative> is $f(y_i\mid\theta_i,\phi)[y_i-b'(\theta_i)]/\phi$. Consequently
$$
0=\frac1\phi\bigl(EY_i-b'(\theta_i)\bigr),\qquad\boxed{\mu_i=b'(\theta_i).}
$$
Differentiating the expectation once more gives $d\mu_i/d\theta_i=\operatorname{Var}(Y_i)/\phi$, so
$$
\boxed{V_i=\operatorname{Var}(Y_i)=\phi b''(\theta_i).}
$$
The <log-likelihood> <derivative> with respect to $\theta_i$ is $(y_i-\mu_i)/\phi$. For the <generalized linear model> with linear predictor $\eta_i=\beta^Tx_i=g(\mu_i)$, the <chain rule> gives
$$
\frac{d\theta_i}{d\mu_i}=\frac1{b''(\theta_i)},\qquad\frac{\partial\mu_i}{\partial\beta}=\frac{x_i}{g'(\mu_i)}.
$$
Using independence to add the <log-likelihood> contributions proves the <score function>
$$
\boxed{\frac{\partial\ell}{\partial\beta}=\sum_{i=1}^n\frac{(y_i-\mu_i)x_i}{g'(\mu_i)V_i}.}
$$
Here $V_i$ includes the dispersion parameter $\phi$. Defining instead a unit-dispersion <variance function> $b''(\theta_i)$ would require retaining an explicit $\phi$ in the denominator.