Can I use Leibniz' rule to differentiate when integral is $\int_0^\gamma f(\gamma,s)dF(s)$, instead of $\int_0^\gamma f(\gamma,s)ds$ ($F$ is a CDF) Specifically, suppose I have a CDF $F$ and an integral
$$
\int_0^\gamma (\gamma-s)dF(s)
$$
and I want to find
$$
\frac{d}{d\gamma} \int_0^\gamma (\gamma-s)dF(s)
$$
(Note that the upper limit $\gamma$ indicates that $s=\gamma$)

If the integral was $\int_0^\gamma (\gamma-s)ds$ then I can apply Liebniz rule. Can I still use Liebniz rule even though my integral is of the form $\int_0^\gamma f(\gamma,s)dF(s)$, instead of $\int_0^\gamma f(\gamma,s)ds$

(if $F$ is differentiable than I can (rigorously/correctly) write $dF(s) = f(s)ds$ and have the standard form. But if $F$ has mass point(s)... then I don't know)
 A: I will interpret the integral as Riemann-Stieltjes integral. Then by integration by parts,
\begin{align*}
\int_{0}^{\gamma} (\gamma-s) \, dF(s)
&= \left[ (\gamma-s)F(s) \right]_{0}^{\gamma} + \int_{0}^{\gamma} F(s) \, ds \\
&= -\gamma F(0) + \int_{0}^{\gamma} F(s) \, ds \\
&= \int_{0}^{\gamma} (F(s) - F(0)) \, ds.
\end{align*}
So we have
$$ \frac{d}{d\gamma} \int_{0}^{\gamma} (\gamma-s) \, dF(s) = F(\gamma) - F(0) = \int_{0}^{\gamma} dF(s) $$
at every continuity point $\gamma$ of $F$. But also notice that this is exactly what we expect when applying the Leibniz integral rule:
$$ \frac{d}{d\gamma} \int_{0}^{\gamma} (\gamma-s) \, dF(s)
\quad``\,=\text{''}\quad \underbrace{(\gamma - \gamma)F'(\gamma)}_{=0} + \int_{0}^{\gamma} dF(s).  $$
A: I am not an expert in probability theory, but I think my knowledge is just sufficient to tackle this problem.
Assume $\gamma>0$ for convenience.
Let $\{a_i\}$ where $a_i>0$ be the discontinuities of $F(s)$. For $a_k<\gamma<a_{k+1}$, we can rewrite your integral as
$$\int_0^{a_1-\epsilon}(\gamma-s)f(s)ds+\left(\sum_{j=1}^{k-1}\int^{a_{j+1}-\epsilon}_{a_j+\epsilon}(\gamma-s)f(s)ds\right)+\int^\gamma_{a_k+\epsilon}(\gamma-s)f(s)ds$$
Differentiating it yields
$$\int_0^{a_1-\epsilon}f(s)ds+\left(\sum_{j=1}^{k-1}\int^{a_{j+1}-\epsilon}_{a_j+\epsilon}f(s)ds\right)+\int^\gamma_{a_k+\epsilon}f(s)ds$$
Although a $\gamma$ appears in the integration limit of the last integral, but if you apply Leibniz integral rule carefully, you can see directly bringing the differentiation into the integral would give the correct result.
EDIT:
I should have explicitly state that $\epsilon$ is to be taken the limit $\to 0^+$.
Also, it might be counterintuitive that the value of the integrand at a point does not affect the value of the integral.
Proof:
Indeed, when decomposing your integral into partitions, I didn’t write out some integrals, whose absolute value is
$$\vert\int^{a_j+\epsilon}_{a_j-\epsilon}f(s)ds\vert\le(a_j+\epsilon-(a_j-\epsilon))|\text{max}_{[a_j-\epsilon, a_j+\epsilon]}f(s)|=2\epsilon |\text{max}_{[a_j-\epsilon, a_j+\epsilon]}f(s)|$$ which vanishes in the limit, as long as the function is not infinite in the region.
