Plot profile likelihood So, I have a generated data from a nonlinear model $$y_i = \eta_1 - 2\theta \eta_2 x_i + \eta_2 x_i^2 + e_i,$$ where $e_i \sim N(0,\sigma^2)$. What I want is to find the profile log-likelihood for $\theta$. Now, like in linear regression, I can write the log-likelihood
$$L(\eta_1,\theta,\eta_2) = -n\log(\sigma) - \frac{n}{2}\log(2\pi) - \frac{1}{2\sigma^2}\sum_{i=1}^{n} \left(y_i - (\eta_1 - 2\theta \eta_2 x_i + \eta_2 x_i^2) \right)^2.$$ Hence, what I need is $\hat{\eta}_1,\hat{\eta}_2$ by maximizing $L(\eta_1,\theta,\eta_2)$ with $\theta$ fixed. But, can I do this analytically? Those maximums should be a function of $\theta$ and also the data (in order to plot the profile log-likelihood for $\theta$)? Thanks
 A: I'll use $\mu_i = \eta_1 - 2\theta\eta_2x_i + \eta_2 x_i^2$ for convenience. If we're thinking of $\mu_i$ as a function of $\theta$, so only $\eta_1$ and $\eta_2$ are parameters, then we can write this as
$$
\mu_i = \eta_1 + \eta_2(-2\theta x_i + x_i^2) = \eta_1 + \eta_2 z_i
$$
for $z_i = -2\theta x_i + x_i^2$. This is just a simple linear regression now so
$$
\begin{aligned}
&\hat\eta_1 = \bar y - \hat\eta_2 \bar z \\&
\hat\eta_2 = \frac{\sum_i (z_i  -\bar z)(y_i - \bar y)}{\sum_i (z_i - \bar z)^2}
\end{aligned}
$$
so all together the profiled log-likelihood is
$$
\ell_p(\theta) = \ell(\hat\eta_1(\theta), \theta, \hat\eta_2(\theta)) \\
= -\frac n2 \log 2\pi\sigma^2 - \frac 1{2\sigma^2}\sum_{i=1}^n (y_i - \hat \eta_1(\theta) - \hat \eta_2(\theta)\cdot(- 2\theta x_i + x_i^2))^2.
$$
Here's an example in R:
set.seed(132)
theta <- 1.23; eta1 <- -.55; eta2 <- .761
sigma <- .234
n <- 500
x <- rnorm(n, -.5)
y <- eta1 - 2 * theta * eta2 * x + eta2 * x^2 + rnorm(n, 0, sigma)

profloglik <- function(theta, sigma, x, y) {
  z <- -2 * theta * x + x^2  # creating the new feature in terms of theta
  mod <- lm(y ~ z)  # using `lm` to do the simple linear regression
  sum(dnorm(y, fitted(mod), sigma, log=TRUE)) # log likelihood 
}

theta_seq <- seq(-10, 10, length=500)
liks <- sapply(theta_seq, profloglik, sigma=sigma, x=x, y=y)

plot(liks ~ theta_seq, type="l", lwd=2,
     main=bquote("Profiled log-likelihood for" ~ theta),
     ylab="profiled log lik", xlab=bquote(theta))


