Double factorial as a sum

I believe the following equality to hold for all integer $l\geq 1$ $$(2l+1)!2^l\sum_{k=0}^l\frac{(-1)^k(l-k)!}{k!(2l-2k+1)!4^k}=(-1)^l(2l-1)!!$$ (it's correct for at least $l=1,2,3,4$), but cannot prove it. I've tried induction but with no success. Any ideas are very welcome.

EDIT: I have come up with an indirect proof. I am still interested in a direct proof. Let $H_m$ denote the $m$th Hermite polynomial defined as (see here for the definition and properties of Hermite polynomials) \begin{align*} H_m(x):=(-1)^m\mathrm{e}^{\frac{x^2}{2}}\frac{d^m}{dx^m}\mathrm{e}^{-\frac{x^2}{2}}, \quad x\in\mathbb{R}. \end{align*}

Let $X$ be a standard Gaussian random variable and $\mathrm{sgn}$ denote the sign function. One can evaluate for all $p\geq 0$ \begin{align} \mathbf{E}\mathrm{sgn}(X)X^{p}&=\left\{\begin{array}{ll} \mathbf{E}|X|^p & \text{if }p \text{ is odd},\\ 0 & \text{if } p \text{ is even}, \end{array}\right. =\left\{\begin{array}{ll} \frac{2^{\frac{p}{2}}\Gamma\left(\frac{p+1}{2}\right)} {\sqrt{\pi}} & \text{if }p \text{ is odd},\\ 0 & \text{if } p \text{ is even}. \end{array}\right.(1) \end{align}

We have the following explicit formula for the Hermite polynomials: for all $m$ and $x\in\mathbb{R}$ \begin{align*} H_m(x)=m!\sum_{k=0}^{\lfloor m/2\rfloor}\frac{(-1)^k}{k!(m-2k)!2^k} x^{m-2k}. \end{align*}

By $(1)$ for all $m=2l+1$ with some $l\geq 1$ we can calculate \begin{align*} a_m&:=\frac{1}{2}\mathbf{E} \mathrm{sgn}(X)H_m(X)= \frac{(2l+1)!}{2\sqrt{\pi}}\sum_{k=0}^l \frac{(-1)^k}{k!(2l-2k+1)!2^k}2^{l-k+1/2}\Gamma(l-k+1)\\ &=\frac{(2l+1)!2^l}{\sqrt{2\pi}}\sum_{k=0}^l \frac{(-1)^k(l-k)!}{k!(2l-2k+1)!4^k}. \end{align*} On the other hand for any $\varepsilon>0$ \begin{align} \sqrt{2\pi}a_m&=\int_0^{\infty}\mathrm{sgn}(x)H_m(x) \mathrm{e}^{-\frac{x^2}{2}}dx\\ &=\int_0^{\varepsilon}\mathrm{sgn}(x)H_m(x) \mathrm{e}^{-\frac{x^2}{2}}dx+(-1)^m\int_{\varepsilon}^{\infty}\frac{d^m}{dx^m} \mathrm{e}^{-\frac{x^2}{2}}dx\notag\\ &=\int_0^{\varepsilon}\mathrm{sgn}(x)H_m(x) \mathrm{e}^{-\frac{x^2}{2}}dx+(-1)^m\left(\frac{d^{m-1}}{dx^{m-1}} \mathrm{e}^{-\frac{x^2}{2}}\right)\bigg|^{\infty}_{\varepsilon}.\quad(2) \end{align} Now \begin{align*} \left(\frac{d^{m-1}}{dx^{m-1}} \mathrm{e}^{-\frac{x^2}{2}}\right)\bigg|^{\infty}_{\varepsilon}= (-1)^{m-1}\mathrm{e}^{-\frac{x^2}{2}}H_{m-1}(x)\bigg|^{\infty}_{\varepsilon}= (-1)^m\mathrm{e}^{-\frac{\varepsilon^2}{2}}H_{m-1}(\varepsilon). \end{align*} Since the integrand in $(2)$ is bounded on $[0,\varepsilon]$ (e.g. here) and $\varepsilon$ is arbitrary we conclude that \begin{align*} \sqrt{2\pi}a_m=H_{m-1}(0)=(-1)^{(m-1)/2}(m-2)!! \end{align*} and thus the desired identity holds for all $l\geq 1$.

We show the identity \begin{align*} \color{blue}{(2l+1)!2^l\sum_{k=0}^l\frac{(-1)^k(l-k)!}{k!(2l-2k+1)!4^k}=(-1)^l(2l-1)!!\qquad\qquad l\geq 1}\tag{1} \end{align*} is valid. In order to do so we divide (1) by the RHS and show the resulting expression is equal to $1$.

We obtain \begin{align*} &\color{blue}{\frac{(2l+1)!2^l(-1)^l}{(2l-1)!!}}\color{blue}{\sum_{k=0}^l\frac{(-1)^k(l-k)!}{k!(2l-2k+1)!4^k}}\\ &\qquad\quad=\frac{(2l+1)!2^l(2l)!!}{(2l)!}\sum_{k=0}^l\frac{(-1)^kk!}{(l-k)!(2k+1)!4^{l-k}}\tag{2.1} \\ &\qquad\quad=(2l+1)\sum_{k=0}^l\binom{l}{k}\binom{2k}{k}^{-1}\frac{(-4)^k}{2k+1}\tag{2.2}\\ &\qquad\quad=(2l+1)\sum_{k=0}^l\binom{l}{k}\int_{0}^1z^k(1-z)^k\,dz(-4)^k\tag{2.3}\\ &\qquad\quad=(2l+1)\int_{0}^1\sum_{k=0}^l\binom{l}{k}(-4z)^k(1-z)^k\,dz\\ &\qquad\quad=(2l+1)\int_{0}^1\left(1-4z(1-z)\right)^k\,dz\tag{2.4}\\ &\qquad\quad=(2l+1)\int_{0}^1(1-2z)^{2l}\,dz\\ &\qquad\quad=\left.(1-2z)^{2l+1}\left(-\frac{1}{2}\right)\right|_{0}^1\tag{2.5}\\ &\qquad\quad=(-1)^{2l+1}\left(-\frac{1}{2}\right)+\frac{1}{2}\\ &\qquad\quad\,\,\color{blue}{=1} \end{align*} and the claim (1) follows.

Comment:

• In (2.1) we use the identity $(2l)!=(2l)!!(2l-1)!!$ and change the order of summation $k\to l-k$.

• In (2.2) we apply $(2l)!!=2^ll!$, do some cancellation and write the expression using binomial coefficients.

• In (2.3) we write the reciprocal of a binomial coefficient using the beta function \begin{align*} \binom{n}{k}^{-1}=(n+1)\int_0^1z^k(1-z)^{n-k}\,dz \end{align*}

• In (2.4) we apply the binomial theorem.

• In (2.5) we do the integration.