How do I prove that the Riemann curvature tensor is indeed a tensor? I know that to prove that something is a tensor I have to show that this thing transforms like a tensor, i.e., like below:
$$ R^{'\alpha}_{\gamma\phi\lambda} = \partial_{\beta}x^{\alpha'}\partial_{\sigma'}x^{\gamma}\partial_{\mu'}x^{\phi}\partial_{\nu'}x^{\lambda}R^{\beta}_{\gamma\phi\lambda}.$$
But I don't know actually how to start the proof. Do I need to apply this tensor to a vector or something like that? Thanks in advance.
 A: The way you prove that something is a tensor really depends on how you define that "something". So, in order to start your proof, you need to bring up the definitions in the first place.
Unfortunately, even the identity that you provide as a criterion of being a tensor is wrong: the indexes must match on the both sides, at least.
Let me change the notation slightly in order to avoid using primes (like $x'$), and instead of your primed quantities I will use just another set of letters.
What really this transformation identity should tell us, is how a quantity transforms when the coordinate system changes from $(x^i)_{i=1,\dots,n}$ to $(y^{\alpha})_{\alpha=1,\dots,n}$ (so, we assume that the dimension of the underlying manifold is $n$). I will refer to these coordinate systems as $(x)$, and $(y)$ respectively, for the sake brevity.
The change of coordinates from $(x)$ to $(y)$ results in a linear transformation of the coordinate local frames from $(\partial_i)_{i=1,\dots,n}$ to $(\partial_{\alpha})_{\alpha=1,\dots,n}$, which again will be referred to as $(\partial_i)$ to $(\partial_{\alpha})$ to save space. This linear transformation (the Jacobian) is represented by a $(n \times n)$-matrix, which in frames $(\partial_i)$ to $(\partial_{\alpha})$ is denoted as $(\partial_i y^{\alpha})$. Since coordinate changes correspond to some diffeomorphisms, the Jacobian in this case is invertible, and its inverse is denoted by $(\partial_{\alpha} y^i)$. Here $i=1,\dots,n$ and $\alpha=1,\dots,n$.
The quantity $R$ under consideration is represented by a set of components $R_{i j}{}^{k}{}_{l}$ in a given coordinate system $(x^i)_{i=1,\dots,n}$. If at any given point $R$ were a multilinear operator, acting on (co-)tangent spaces, then its components in the aforementioned coordinates frames (recall, that a frame is just a basis of a vector space, smoothly varying along the manifold), would be given by a formula like
$$
R_{\alpha \beta}{}^{\gamma}{}_{\delta} = \partial_{\alpha} y^i \, \partial_{\beta} y^j \, \partial_k y^{\gamma} \, \partial_{\delta} y^l \, R_{i j}{}^{k}{}_l \tag{$\star$}
$$
where for each repeated pair of indexes we assume the Einstein summation. This is pure multilinear algebra. Conversely, if the identity $(\star)$ holds, then $R$ is clearly a multilinear operator.
Now, if you had an explicit formula for $R$, you could probably just read off from that formula, whether $R$ is multilinear or not.
In an attempt to answer another question, I have presented an explicit formula for the second covariant derivative of a vector in terms of the Christoffel symbols:
$$
\begin{align}
(\nabla \nabla X)_{a b}{}^{c} & = \partial_a \partial_b X^c - \Gamma_{a b}{}^{d} \Gamma_{d e}{}^{c} X^e +  \Gamma_{a d}{}^{c} \Gamma_{b e}{}^{d} X^e  \\
& + (\partial_a \Gamma_{b d}{}^{c}) X^d + \Gamma_{b d}{}^{c} \partial_a X^d - \Gamma_{a b}{}^{d}(\partial_d X^c)  X^e + \Gamma_{a d}{}^{c} (\partial_b X^d) 
\end{align}
$$
Examining the above formula, we can notice that anti-symmetryzing the indices $a$ and $b$, the terms with the partial derivative of the components of vector $X$ disappear, and we are left with the expression
$$
2(\nabla \nabla X)_{[a b]}{}^{c} = 2 (\partial_{[a} \Gamma_{b] d}{}^{c}) X^d + 2 \Gamma_{[a | d |}{}^{c} \Gamma_{b] e}{}^{d} X^e  \\
$$
that, after renaming indexes in the second term, we can rewrite as
$$
2(\nabla \nabla X)_{[a b]}{}^{c} = \big( 2 (\partial_{[a} \Gamma_{b] d}{}^{c}) + 2 \Gamma_{[a | f |}{}^{c} \Gamma_{b] d}{}^{f} \big) X^d  \\
$$
In the RHS we see a linear combination of partial derivative of the metric tensor, and this combination acts linearly onto the components of vector $X$. Therefore, this combination is in fact a set of components of some multilinear operator, which is denoted as 
$$
R_{a b}{}^{c}{}_{d} := \partial_{a} \Gamma_{b d}{}^{c} - \partial_{b} \Gamma_{a d}{}^{c} + \Gamma_{a e}{}^{c} \Gamma_{b d}{}^{e} - \Gamma_{b e}{}^{c} \Gamma_{a d}{}^{e}
$$
and referred to as the Riemann curvature operator.
