Write down the log-likelihood function for this model 
Consider the regression model
$$ Y_i = ax_i^3 + \epsilon_i, \hspace{1cm} i = 1,...,n$$
with $\epsilon_i \sim N(0, \sigma^2)$ and $\epsilon_i, \epsilon_j$ independent for $i \neq j$. Write down the log likelihood function for this model.

$a$ is another variable, not a constant. For some reason in the answers, they use the PDF for the normal distribution to do this. Why? How do you know that's what you use and not just the regression model?
EDIT: This is revision, a question in a past paper. (That is, by the tag description, this is homework.)
 A: The only way a maximum likelihood estimation can be implemented here (as @whuber's comments make clear) is to assume that $a$ is a "variable" in the sense of being a "varying coefficient" or a "random coefficient". Then there are various ways to model this situation (the fact that $x_i$ is raised to the third power does not impede us). For example, an (enriched, and for a reason, see at the end) variant of the Hildreth-Houck (1968) model of random coefficients is  
$$y_i = a_iz_i + \varepsilon_i$$,
with
$$ a_i = a + u_i,\; u_i\sim N(0,\sigma^2_u),\;\; \varepsilon_i \sim N(0,\sigma^2_{\varepsilon})$$
with the two errors independent, and $z_i$ exogenous to both. So the model can be re-written
$$y_i = az_i + w_i,\qquad w_i = z_iu_i + \varepsilon_i,\qquad w_i \sim N(0, s^2_i),\;\; s^2_i = z_i^2\sigma^2_u+\sigma^2_{\varepsilon}$$
The model can be seen as a standard linear regression model with heteroskedastic errors.
The log-likelihood function here is ($\mathbf q = \{a,\sigma^2_u, \sigma^2_{\varepsilon}\}$)
$$\ell (\mathbf q \mid \mathbf z) = k -\sum_{i=1}^n\ln s^2_i -\frac 12\sum_{i=1}^n\left(\frac {(y_i-az_i)^2}{s^2_i}\right)$$
To find the MLE for $a$ we have
$$\frac {\partial \ell (\mathbf q \mid \mathbf z)}{\partial a} = \sum_{i=1}^n\left(\frac {(y_i-az_i)z_i}{s^2_i}\right) = \sum_{i=1}^n\left(\frac {y_iz_i}{s^2_i}\right)-a\sum_{i=1}^n\left(\frac {z_i^2}{s^2_i}\right)= 0$$
$$\Rightarrow \hat a_{MLE} = \left[\sum_{i=1}^n\left(\frac {z_i^2}{s^2_i}\right)\right]^{-1}\cdot\sum_{i=1}^n\left(\frac {y_iz_i}{s^2_i}\right)$$
Now, since $y_iz_i = (az_i + w_i)z_i = (az_i + z_iu_i + \varepsilon_i)z_i = az_i^2 + (z_i^2u_i + \varepsilon_iz_i)$ we obtain
$$\hat a_{MLE} = a + \left[\sum_{i=1}^n\left(\frac {z_i^2}{z_i^2\sigma^2_u+\sigma^2_{\varepsilon}}\right)\right]^{-1}\cdot\sum_{i=1}^n\left(\frac {z_i^2u_i + \varepsilon_iz_i}{z_i^2\sigma^2_u+\sigma^2_{\varepsilon}}\right)$$
Under the exogeneity assumptions, the MLE is consistent.
We close by the following note: the original Hildreth-Houck model did not assume the existence of two random disturbances, i.e it assumed $\sigma^2_{\varepsilon} =0$ (assuming two errors was the enrichment mentioned previously). Then in order for the MLE to be consistent, we need the regressor to be bounded away from zero for all observations, because otherwise for some $i$ the variance of $w_i$ which would now be $s^2_i = z_i^2\sigma^2_u$ would equal zero for this $i$ and the value of the likelihood would explode. See, for example Zaman (2002)
A: In terms of probability theory, this is a application of a technique known as change of variables.  This tells you how the PDF for one random variable (in your case $\epsilon_i$) is related the PDF for a function of that random variable ($Y_i$ in your case).  The general rule is given as
$$f(Y_i)=f(\epsilon_i)|\frac{\partial\epsilon_i}{\partial Y_i}|$$
This is slightly slack notation, but this is usually how the result is presented.  This is also equivalent to the change of variables technique used for integrals $\int f(\epsilon_i)d\epsilon_i$.  Note that for the linear regression model the derivative term is one.  This means the PDF for $Y_i$ is given by a simple substitution into the PDF for $\epsilon_i$ via the (rearranged) model equation $\epsilon_i=Y_i-ax_i^3$.  Upon making this substitution we find that the PDF for $Y_i$ is a normal distribution with mean $ax_i^3$ and variance $\sigma^2$.
