As far as I can tell, you have two questions. I will answer them separately. Throughout I will set $\hbar = 1$.
1) Where does the sum $\sum_{n = -\infty}^{\infty}$ come from?
The standard path integral form for the propagator of a free particle is
$$ \langle x_{2} | e^{-\beta H} | x_{1} \rangle = \int \mathcal{D} x(\tau) \exp(-S_{E}[x(\tau)]),$$
where the integrals runs over all paths $x(\tau)$ satisfying the boundary conditions $x(0) = x_{1}$ and $x(\beta)=x_{2}$, and $S_{E}$ is the Euclidean action.
The crux of this answer relies on the fact that in the path integral, we integrate over all possible paths satisfying the given boundary conditions. What the change of variables from $x(\tau)$ to $s(\tau)$ does is to make how we count these paths more explicit.
In the question we define a variable $s(\tau)$ that satisfies $s(0) = s(\beta) = 0$. Other than these boundary conditions, $s(\tau)$ is completely unrestricted (i.e. it is not restricted to a circular topology). Now, for each choice of $s(\tau)$, we can define a path
$$ x_{n}(\tau) = x_{1} + \dfrac{\tau}{\beta} \left[(x_{2}-x_{1}) + n L \right] + s(\tau), $$
for some integer $n$, and this will be a valid path for the path integral, since the circular topology means that we identify points modulo $L$. Furthermore, between the choices for $s(\tau)$ and $n$, we have enumerated all possible paths that satisfy the given boundary conditions.
Thus, in the expression for the propagator, we have to integrate over $s(\tau)$ and sum over $n$ to make sure we account for all possible paths $x(\tau)$. This gives us the expression
$$ \langle x_{2} | e^{-\beta H} | x_{1} \rangle = \sum_{n} \int \mathcal{D} s(\tau) \exp(-S_{E}[x_{n}(\tau)]),$$
with $x_{n}(\tau)$ given as above.
2) Why should $Z_{0}$ be the diagonal matrix element of the propagator on $\mathbb{R}$?
As you stated in your question, the Euclidean action $S_{E}[x_{n}(\tau)]$ separates into two parts, given by
$$ S_{E}[x_{n}(\tau)] = S_{n} + S_{E}[s(\tau)], $$
where $S_{n} = \frac{m}{2 \beta} \left[ (x_{2} - x_{1}) + n L \right]^{2}$. This separation is helpful because it separates the dependence on $n$ and $s(\tau)$.
Since $S_{n}$ is a constant, we can factor it out of the path integral. This part gives the summand in the expression for the propagator that you are trying to derive.
The remaining part is
$$ \int \mathcal{D} s(\tau) \exp(-S_{E}[s(\tau)]),$$
where we recall that the path $s(\tau)$ has the boundary conditions $s(0) = s(\beta) = 0$, but is otherwise unrestricted (i.e. it is not restricted to a circular topology). But this is precisely what we'd write for the propagator of a free particle in one dimension moving from $x=0$ to $x=0$ in imaginary time $\beta$! And from translational invariance, this shouldn't depend on the fact that the particle starts and finishes at $x=0$, just the fact that the path is periodic. Hence the remaining part gives us $Z_{0} = \langle x | e^{-\beta H} | x \rangle$ for a free particle moving on the real line. This is the "diagonal" matrix element of the propagator because it is simply the propagator $ \langle x_{2} | e^{-\beta H} | x_{1} \rangle $ evaluated at $x_{1} = x_{2} = x$.
