Nonlinear regression: What's the best estimator for this problem? Context and goals:
Consider $x$ a scalar (deterministic) independent variable and $y = f(x;\beta)+\eta$ be a dependent random variable obtained trough $f$ with some parameters $\beta$ and $\eta$ is some noise modeled as a Gaussian with zero mean and unknown variance. Concretely, chose
$$
f(x;\beta) = \beta_1 + \frac{\beta_2}{1+e^{\beta_3x}}
$$
Now, assume that we made $M$ experiments (with $M$ presumably low) and obtained $M$ pairs of measurements $(x_1,y_1),\dots,(x_M,y_M)$.

*

*Question: Given the data, How to obtain a way to estimate of $y$ for any $x$ as well as its "confidence"?

Let me explain my attempt to solve this problem:
My idea was to fit the parameters to the data using nonlinear regression by means of optimization (gradient descent) for example. Thus obtain an approximate model $f(\bullet;\hat{\beta})$. Then, for certain value of $x$ we can model that its corresponding $y$ is a Gaussian centered at $f(\bullet;\hat{\beta})$ (right?). However, since I don't know the variance of the noise $\eta$, the most we can do is to compute the sample variance of $y_1,\dots,y_M$ and use a student's-t distribution instead.
The problem is that, I'm not totally convinced that $f(\bullet;\hat{\beta})$ can be used as the mean of $y$. I guess one has to check if $f(\bullet;\hat{\beta})$ is an unbiased estimate of $y$. However, I don't know how to check that. Moreover, by using the sample variance I only get a number which represents the "confidence" for all $y$ which seems unsatisfactory to me since I would expect to be more confident of values of $y$ for $x$ near the samples $x_1,\dots,x_M$. But maybe my intuition is not correct.
 A: Given the nonlinear regression model
$$
Y_{i} = f( x_{i}, \underline{\beta} ) + \eta_{i} \quad \textrm{where} \quad  \eta_i \overset{\text{iid}}{\sim} \mathcal{N}(0, \sigma^{2})
$$
and given the data
$$
(x_{i}, y_{i}), i = 1, ..., M
$$
a standard approach to estimate the unknown parameters $\underline{\beta} = (\beta_{1}, ..., \beta_{p})$, is to apply the method of least squares by minimizing the sum of squared residuals
$$
\sum_{i=1}^{M} \left( y_{i} - f(x_{i}, \underline{\beta}) \right)^{2}
$$
Since this problem is not linear in the unknowns $\beta_{1}, ..., \beta_{p}$, the estimation can not be derived explicitly as in linear regression, but we can solve it with the help of an iterative numeric procedure such as Gauss-Newton method, see: http://en.wikipedia.org/wiki/Gauss-Newton_method.
The least squares estimator $\hat{\beta}$ is asymptotically multivariate normally distributed with mean $\beta$ and therefore unbiased.
For a given $x_{0}$, the function $f$ can be approximated with:
$$
f(x_{0}, \underline{\hat{\beta}}) \approx f(x_{0}, \underline{\beta})
+ \frac{\partial f(x_0, \underline{\beta}) }{ \partial \underline{\beta}}
(\underline{\hat{\beta}} - \underline{\beta})
$$
and the approximate $(1 - \alpha)$ confidence interval for the function value $f(x_{0}, \underline{\beta})$ is then
$$
f(x_{0}, \underline{\hat{\beta}}) \pm q_{1-\alpha/2}^{t_{m - p}} \cdot 
\textrm{se}( f(x_{0}, \underline{\hat{\beta}}) )
$$
We can estimate the expected value $E(Y_{0}) = f(x_{0}, \underline{\beta})$ with $f(x_{0}, \underline{\hat{\beta}})$. And to get the interval where the observations will lie with high probability we also have to take into account the variance of the error $\eta_{0}$
$$
f(x_{0}, \underline{\hat{\beta}}) \pm q_{1-\alpha/2}^{t_{m - p}} \cdot 
\sqrt{
\hat{\sigma}^{2} + 
\textrm{se}( f(x_{0}, \underline{\hat{\beta}}) )^{2}
}
$$
This is called prediction interval or forecast interval.
