Implements an iterative mean-variance panel regression estimator
that allows both the mean and variance of the dependent variable to be
functions of covariates. The method alternates between estimating a mean
equation (using generalized linear models with Gaussian family) and a
variance equation (using generalized linear models with Gamma family on
squared within-group residuals) until convergence. Based on the methodology
in Mooi-Reci and Liao (2025)