Is a probability density function necessarily a $L^2$ function? If a nonnegative continuous real valued function $f$ is integrable over $\mathbb{R}$ with
$$\int_\mathbb{R} f\,\mathrm{d}x = 1,$$ 
does it hold true
$$\int_\mathbb{R} f^2 \,\mathrm{d}x<\infty?$$
Motivation: I am wondering if the mean squared error (MSE) is well defined, since we need the target density function to be in $L^2$ space in order to have a finite MSE.
 A: The answer is no. But your continuity requirement makes it more difficult to find an explicit (constructed) example.
Define the function $\varphi$ on $\mathbb{R}$ as:
$$\forall x\in\mathbb{R},\quad\varphi(x)=\begin{cases}\sin(x)&\text{if $x\in[0,\pi]$}\\0&\text{otherwise.}\end{cases}$$
Clearly, the function $\varphi$ is continuous on $\mathbb{R}$ (it's just a single arch of sine). It is well-known that $\lVert\varphi\rVert_1=2$ and that $\lVert\varphi\rVert_2^2=\pi/2$.
Define the sequence of function $(\varphi_n)_{n\geq1}$ on $\mathbb{R}$ by:
$$\forall n\geq1,\quad\forall x\in\mathbb{R},\quad\varphi_n(x)=\varphi\bigl(8^n(x-2n\pi)\bigr).$$
It should be clear that the support of the $\varphi_n$'s are all distinct: the support of $\varphi_n$ is $(2n\pi,8^{-n}\pi+2n\pi)$.
Moreover, for all $n\geq1$,
$$
\lVert\varphi_n\rVert_1=\int_{\mathbb{R}}\varphi_n(x)\,\mathrm{d}x=\frac2{8^n}\quad\text{and}\quad\lVert\varphi_n\rVert_2^2=\int_{\mathbb{R}}\varphi_n(x)^2\,\mathrm{d}x=8^{-n}\frac{\pi}2.$$
Now, define the series of functions $f$ by:
$$\forall x\in\mathbb{R},\quad f(x)=\frac12\sum_{n=1}^{+\infty}4^n\varphi_n(x)$$
(the graph of $f$ consists of bumps that are narrower but higher).
It should be clear that $f$ is well-defined and continuous and non-negative.
Since the support of the $\varphi_n$'s are distinct, we can integrate $f$ term by term:
$$\lVert f\rVert_1=\int_{\mathbb{R}}f(x)\,\mathrm{d}x=\frac12\sum_{n=1}^{+\infty}4^n\int_{\mathbb{R}}\varphi_n(x)\,\mathrm{d}x=\sum_{n=1}^{+\infty}\frac{4^n}{8^n}=1.$$
Also, since the support of the $\varphi_n$'s are distinct, for all $x\in\mathbb{R}$:
$$f(x)^2=\frac14\sum_{n=1}^{+\infty}16^n\varphi_n(x)^2$$
and
$$\lVert f\rVert_2^2=\int_{\mathbb{R}}f(x)^2\,\mathrm{d}x
   =\frac14\sum_{n=1}^{+\infty}16^n\lVert\varphi_n\rVert_2^2=\frac14\sum_{n=1}^{+\infty}2^n\frac\pi2=+\infty.$$
