infeasible when the number of site locations increases [@davison+gholamrezaee:2011].

In the remainder of this section, we show that the gap between improved estimation methods and applied modelling using general $M$-estimators can be bridged either through simple empirical bias-reducing adjustments to the estimating functions (\[estimating\_eqautions\]) that do not involve expectations, or empirical bias-reducing penalties to the objective function (\[objective\_function\]).

Family of bias-reducing adjustments to estimating functions
-----------------------------------------------------------

Suppose that the assumptions in Section \[assumptions\] hold and consider the estimator $\tilde{{}\theta}$ that results from the solution of the adjusted estimating equations $$\label{adjusted_estimating_equations}
  \sum_{i=1}^k {{}\psi}^i({{}\theta}) + {{}A}({{}\theta}) = {{}0}_p \, ,$$ where both ${{}A}({{}\theta}) = {{}A}({{}\theta}, {{}Y}, X)$ and its derivatives with respect to ${{}\theta}$ are $O_p(1)$ as $n$ grows.

Using assumptions \[consistency\]-\[derivative\_orders\] and index notation, with the indices taking values in the set $\{1,\ldots,p\}$, a calculation similar to that in @mccullagh:2018 [Section 7.3] can be used to show that the expansion of ${{}0}_p = \sum_{i=1}^k {{}\psi}^i(\tilde{{}\theta}) + {{}A}(\tilde{{}\theta})$ about $\bar{{}\theta}$ results in a stochastic Taylor expansion for $\tilde{{}\theta} - \bar{{}\theta}$ of the form $$\begin{aligned}
  \label{taylor}
  \tilde{\theta}^{r} - \bar{\theta}^{r} & =  H^{r} + H^{a} H_{a}^{r} + \frac{1}{2} H^{a} H^{b} \mu_{ab}^{r} + A^{r} +\\ \notag
                                   &+ H^{a} H_{a}^{b} H_{b}^{r} + \frac{1}{2} H^{a} H^{b} H_{c}^{r} 