Show that the pdf is ${1 \over {(n-1)}!} \sum\limits_{n \leq j \leq \lfloor x\rfloor}(-1)^j \binom{n}{j}(x-j)^{n-1}$
Solution 1:
Most proofs I found seem to be rather messy so I figured I'd try to give a simpler proof.
Let $f_n$ be the distribution of $\sum_{i=1}^{n} X_i$ with $X_i$ all independent and uniformly distributed on $[0,1]$, also let $F_n$ be the cumulative distribution of $f_n$.
Now, as you discovered your formula works for $f_1$ so we can use induction to show that it works for all other $n$. So let's assume that
$$ f_n = {1 \over {(n-1)!}} \sum_{j = \lfloor x\rfloor}^{n} (-1)^j \binom{n}{j}(x-j)^{n-1} \qquad 0 \leq x \leq n $$
from this it is easy to show that
$$ F_n = {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n} (-1)^j \binom{n}{j}(x-j)^{n} \qquad 0 \leq x \leq n. $$
And some basic probability theory shows that
$$ \def\d{\textrm{d}} f_{n+1}(x) = \int_{-\infty}^{\infty} f_{n}(x-t) f_1(t) \, \d t = \int_0^1 f_{n}(x-t) \, \d t = F_n(x) - F_n(x-1) $$
now note that
$$ \begin{align} F_n(x-1) = {1 \over {n!}} \sum_{j = \lfloor x - 1\rfloor}^{n} (-1)^j \binom{n}{j} (x-1-j)^n\\ = - {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n+1} (-1)^{j} \binom{n}{j-1} (x-j)^n \end{align} $$
hence (using the fact that $\binom{n}{n+1} = 0$)
$$ \begin{align} f_{n+1}(x) = F_n(x) - F_n(x-1) &= {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n+1} (-1)^{j} \binom{n}{j} (x-j)^n + {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n+1} (-1)^{j} \binom{n}{j-1} (x-j)^n\\ &= {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n+1} (-1)^{j} \left(\binom{n}{j} + \binom{n}{j-1}\right) (x-j)^n\\ &= {1 \over {n!}} \sum_{j = \lfloor x\rfloor}^{n+1} (-1)^{j} \binom{n+1}{j} (x-j)^n \end{align} $$
finishing the proof.
Solution 2:
Let the random variable be $y=\sum_{i=1}^n x_i$. I will prove this using Laplace transform. A very neat trick is in general if $X_1, \cdots, X_n$ are radom variables with a joint probability distribution $p(x_1, \cdots, x_n)$, $f$ any function of $n$ variables and $y=f(x_1, \cdots, x_n)$ is any new random variable, then $P(y)$ the pdf of $y$ is $P(y)=\left\langle \delta(y-f(x_1, \cdots, x_n)\right\rangle$, where $\delta$ is the $\delta$ function and $\langle \bullet \rangle$ is the expectation value via $p(x_1, \cdots, x_n)$. Therefore in our case $$ P(y)=\left\langle \delta\left(y-\sum_{i=1}^n x_i\right)\right\rangle = \int_0^1 dx_1dx_2\cdots dx_n \delta\left(y-\sum_{i=1}^n x_i\right) $$ Doing a Laplace transform $$ \begin{aligned} \tilde{P}(s) &=\int_{0}^\infty dy \: e^{-sy} \int_0^1 dx_1dx_2\cdots dx_n \delta\left(y-\sum_{i=1}^n x_i\right)= \left(\int_0^1 dx e^{-sx}\right)^n\\&=\left[\frac{1-e^{-s}}{s}\right]^n =\sum_{j=0}^n {n\choose j}(-1)^j \frac{e^{-sj}}{s^n} \end{aligned} $$ At the same time, if $\mathcal{L}$ denotes the Laplace transform and $\theta(x)$ is the step function, and $\alpha$ any non-negative real number $$ \begin{aligned} \mathcal{L}[(y-\alpha)^{n-1}\theta( y-\alpha)]&=\int_0^{\infty} dy (y-\alpha)^{n-1} e^{-sy}\theta(y-\alpha) \\ .\xrightarrow{z=y-\alpha}&\: e^{-s\alpha}\int_{-\alpha}^\infty dz \:z^{n-1} e^{-sz}\:\theta(z)=e^{-s\alpha}\int_{0}^\infty dz \:z^{n-1} e^{-sz}=(n-1)!\:\frac{e^{-s\alpha}}{s^n} \end{aligned} $$ where in the last step we have used the definition of Gamma function. Henceforth $$ \tilde{P}(s)=\mathcal{L}\left[\frac{1}{(n-1)!}\sum_{j=0}^n {n\choose j}(-1)^j (y-j)^{n-1}\theta(y-j)\right] $$ meaning (by inverse Laplace transform) $$ \boxed{ P(y)=\frac{1}{(n-1)!}\sum_{j=0}^{\lfloor y\rfloor} {n\choose j}(-1)^j (y-j)^{n-1}} $$