The short answer is no, your sample estimates will not be zero. But, they will converge to zero as you get more data. Details follow.
Anyone is welcome to refine the convergence claims that I make here, but:
Suppose your model says $Y = X\beta + \epsilon$ with $X$ being a fixed design matrix of full column rank and $\epsilon$ being i.i.d. Gaussian noise. Suppose we have $n$ samples. Then the MLE for $\beta$ is $$\hat \beta = (X^TX)^{-1}X^TY = (X^TX)^{-1}X^TX\beta + (X^TX)^{-1}X^TX\epsilon = \beta + (X^TX)^{-1}X^T\epsilon $$, and it converges to the true value under "reasonable" asymptotic scenarios. You are correct that the estimates need not satisfy the same constraints that govern the true parameters -- to see this, draw an iid standard Normal sample and verify that the mean differs from zero. The closest claims I could justify:
- If $f$ is linear, then $E[f(\hat \beta)]$ is zero in expectation, as you noticed.
- If $f$ is continuous, then $f(\hat \beta)$ converges to zero in probability.
- If $f$ is twice differentiable with gradient of $\nabla f(\beta) = D$, then $f(\hat \beta)$ converges to zero and the distance from zero has a variance of approximately $D^T(X^TX)^{-1}D$ (via the delta method). Footnote: the limiting variance usually has a factor of $1/n$; in my case, that's buried inside the $(X^TX)^{-1}$. For instance, if we are estimating the mean, then $X$ is a column of ones and $(X^TX)^{-1} = 1/n$.
If you really need those constraints, you should try to impose them, or if you just kinda want them, perhaps you could add a penalty such as $-\lambda f^2(\beta)$ to your likelihood.
EDIT: more info on convergence in probability
We are dealing with convergence in probability, so I want to show that for every $q_0, \epsilon > 0$, there exists a sample size such that $P(|f(\hat \beta) - f(\beta)| > \epsilon) < q_0$. My tools of choice: continuity and convergence in probability of $\hat \beta$.
Suppose $\hat \beta$ converges to $\beta$ in probability. That means for every $p_0, \delta > 0$, there exists a sample size such that $P(||\hat \beta - \beta|| > \delta) < p_0$. This is a lot like what I want to show: if I choose $p_0=q_0$, then I just need to choose $\delta$ such that the "bad behavior in $f$" event, $|f(\hat \beta) - f(\beta)| > \epsilon$, is contained in the "bad behavior in $\hat \beta$" event. Fortunately, $f$ is continuous, so for any $\epsilon$ there exists a $\delta$ such that $|f(\hat \beta) - f(\beta)| < \epsilon$ whenever $||\beta - \beta|| < \delta$. Please forgive sloppiness in strict/weak inequalities.
