## Appendix A Theorems on Solving the Score-Matching Equation
\right.\]  
_Furthermore, the solutions are associated with the Lambert \(W\)-function \(P(z)\), which solves \(z=we^{-w}\) (see \(R\) function_ lambert_W in package_ VGAM_[65]). It can be verified that, for cases (ii) and (iii), \(\tilde{\eta}_{k}=P(-\bar{S}e^{\bar{S}})-\bar{S}\); for cases (iv) and (v), \(\tilde{\eta}_{k}=\bar{S}-P(-\bar{S}e^{\bar{S}})\). Nevertheless,_ lambert_ does not provide solutions for cases (i) and (vi) for now. One may use a quasi-Newton algorithm to solve them._  
**Theorem A4**. : _For the Poisson model with a log link, suppose either \(\eta_{i}>0\) for all \(i\in I_{k}\) or \(\eta_{i}<0\) for all \(i\in I_{k}\). Then, there are up to two solutions \(\tilde{\eta}_{k,1}\in[\eta_{k}^{\wedge},\mu(\tilde{y}_{k})]\) and \(\tilde{\eta}_{k,2}\in[u(\tilde{y}_{k}),\eta_{k}^{\vee}]\) solving (3), where \(u(y)\geq-1\) denotes the unique root of the transcendent equation \(e^{\eta}(1+\eta)=y\) given \(y\geq 0\). Special cases include the following:_  
1. _If_ \(y_{i}=0\) _and_ \(\eta_{i}\geq 0\) _for all_ \(i\in I_{k}\)_, then there exists a unique solution_ \(\tilde{\eta}_{k}=1+u(-\bar{S}/e)\geq 0\)_._
2. _If_ \(y_{i}=0\) _and_ \(\eta_{i}\in(-1,0)\) _for all_ \(i\in I_{k}\)_, then there exists a unique solution in_ \((-1,0)\)_._
3. _If_ \(y_{i}=0\) _and_ \(\eta_{i}\leq-1\) _for all_ \(i\in I_{k}\)_, then there exists a unique solution in_ \((-\infty,-1]\)_._  
**Remark A3**. : _Theorem A4 suggests that the observations with \(y_{i}=0\) be extracted and further partitioned into three sub-blocks according to \(\eta_{i}\in(-\infty,-1]\), \((-1,0)\) and \([0,\infty)\), respectively, or with splitting points \(-1\) and \(0\). Then, \(\tilde{\eta}_{k}\) solving (3) within each sub-block is unique and can be obtained easily. Moreover, all three special cases in Theorem A4 are related to the well-known Lambert \(W\)-function \(P(z)\), which solves \(z=we^{zw}\). The function \"lambert_\" from \(R\) package_ VGAM _may be used for this purpose as well._  
**Theorem A5**. : _For the gamma model with a reciprocal link, \(\tilde{\eta}_{k}=n_{k}^{-1}\sum_{i\in I_{k}}\eta_{i}\) is the unique solution solving (3). Furthermore, if \(y_{i}<G(\eta_{i})\) for all \(i\in I_{k}\), then \(\tilde{\eta}_{k}\in(0,\tilde{y}_{k}^{-1})\); if \(y_{i}>G(\eta_{i})\) for all \(i\in I_{k}\), then \(\tilde{\eta}_{k}\in(\tilde{y}_{k}^{-1},\infty)\)._  
**Remark A4**. : _For the gamma model with a reciprocal link, \(\tilde{\eta}_{k}=\tilde{\eta}_{k}=\tilde{\mathbf{X}}_{k}\boldsymbol{\beta}\) is the only solution solving (3), where \(\tilde{\mathbf{X}}_{k}=n_{k}^{-1}\sum_{i\in I_{k}}\mathbf{X}_{i}\) is the mean representative. In this case, \(\
u(\eta)\) is a constant (see Table 1 in [43]), and we have the following:_  
\[\tilde{\mathbf{X}}_{k}=\frac{\sum_{i\in I_{k}}(y_{i}-\eta_{i}^{-1})\mathbf{X} _{i}}{n_{k}(\tilde{y}_{k}-\tilde{\eta}^{-1})}\] (A2)  
_according to (4). That is, \(\tilde{\mathbf{X}}_{k}\
eq\tilde{\mathbf{X}}_{k}\) as in (A2), even if \(\tilde{\eta}_{k}=\tilde{\mathbf{X}}_{k}\boldsymbol{\beta}\)._  
**Theorem A6**. : _For the inverse Gaussian model with an inverse-square link, in general, there are up to two solutions, i.e.,_  
\[\tilde{\eta}_{k,1}=\frac{1-4\tilde{y}_{k}\bar{S}-\sqrt{1-8\tilde{y}_{k}\bar{S}} }{2\tilde{y}_{k}^{2}},\;\;\tilde{\eta}_{k,2}=\frac{1-4\tilde{y}_{k}\bar{S}+ \sqrt{1-8\tilde{y}_{k}\bar{S}}}{2\tilde{y}_{k}^{2}}\]  
_solving (3), where \(\bar{S}=n_{k}^{-1}\sum_{i\in I_{k}}\frac{1}{2}(\eta_{i}^{-1/2}-\tilde{y}_{k}) \eta_{i}\), and \(\tilde{\eta}_{k,1}\leq(4\tilde{y}_{k}^{2})^{-1}\leq\tilde{\eta}_{k,2}\). Furthermore, if \(y_{i}>G(\eta_{i})\) for all \(i\in I_{k}\), then \(\bar{S}<0\) and \(\tilde{\eta}_{k}=(1-4\tilde{y}_{k}\bar{S}+\sqrt{1-8\tilde{y}_{k}\bar{S}})/(2 \tilde{y}_{k}^{2})\) is the only solution._