Assume we have a probablistic model $f_{\theta}(x)$ and try to estimate the parameter $\theta$ based on data $x$ with some procedure that yields a biased estimator $$E[\hat{\theta}]=\theta + \eta,$$ where $\eta$ represents the bias.
As far as I know we can estimate the bias $\eta$ via bootstrap, applying some procedure as follows:
1) Run the model in order to obtain $\hat{\theta}$.
2) Draw $B$ bootstrap samples $X_1,...,X_B$ and re-estimate the model $B$ times in order to obtain paramter estimates $\hat{\theta}_1,...,\hat{\theta}_B$.
3) Estimate the bias via $\hat{\eta}=\hat{\theta}-(n^{-1}\sum_{i=1}^B\hat{\theta}_i-\hat{\theta})$.
However, this procedure implies that the estimate of the bias is $unbiased$ itself. I am, however confronted with a situation where $$E[\hat{\theta}]=\theta + \eta(\theta),$$ and hence the bias a function of the unkown $\theta$ we cannot obtain valid estimates for $\eta$ since simply pluggin-in $\hat{\theta}$ for $\theta$ does not work as we cannot hope to appoximate $\eta(\theta)$ with $\eta(\hat{\theta})$. Do you have any suggestions how to solve such a problem (at least approximately)?