This requires only the inverse of the Cayley transform. Start with
$$
(U-I)=(A-iI)(A+iI)^{-1}-(A+iI)(A+iI)^{-1}=-2i(A+iI)^{-1}.
$$
It follows that $\mathcal{N}(U-I)=\{0\}$ and $\mathcal{R}(U-I)=\mathcal{D}(A)$. Similarly,
$$
(U+I) = 2A(A+iI)^{-1} = iA(U-I).
$$
Let $U=\int_{T}\lambda dF(\lambda)$, and, for each $0 < \delta < \pi$, define $G_{\delta}$ to be the characteristic function of the arc $\{ e^{i\theta} : \theta \in [\delta,2\pi-\delta]\}$. Then
$$
P_{\delta} = \int G_{\delta}(\lambda)dF(\lambda)
$$
is a projection with $P_{\delta}x \in \mathcal{D}(A)$ because
$$
Q_{\delta}=\int_{T} G_{\delta}\frac{1}{\lambda-1}dF(\lambda)
$$
is bounded and $(U-I)Q_{\delta}=P_{\delta}$ implies that the range of $P_{\delta}$ is in $\mathcal{R}(U-I)=\mathcal{D}(A)$. Furthermore,
$$
iAP_{\delta} = iA(U-I)Q_{\delta}=(U+I)Q_{\delta}=\int_{T}G_{\delta}(\lambda)\frac{\lambda+1}{\lambda-1}dF(\lambda) \\
AP_{\delta} = \int_{T}i\frac{1+\lambda}{1-\lambda}G_{\delta}(\lambda)dF(\lambda).
$$
Because $x \in \mathcal{D}(A)$ iff $x = (U-I)y$ for some $y$, then, for all $x \in \mathcal{D}(A)$, one has
$$
\begin{align}
P_{\delta}Ax & = P_{\delta}A(U-I)y \\
& =-iP_{\delta}(U+I)y\\
& =-i(U+I)P_{\delta}y \\
& = A(U-I)P_{\delta}y \\
& = AP_{\delta}(U-I)y = AP_{\delta}x
\end{align}
$$
If $x \in \mathcal{D}(A)$, then
$$
Ax = \lim_{\delta\downarrow 0}P_{\delta}Ax=\lim_{\delta\downarrow 0}AP_{\delta}x
= \lim_{\delta\downarrow 0}\int i\frac{1+\lambda}{1-\lambda}G_{\delta}(\lambda)dF(\lambda)x.
$$
Therefore, by the monotone convergence theorem, if $x\in\mathcal{D}(A)$, then
$$
\begin{align}
\|Ax\|^{2} & =\lim_{\delta\downarrow 0}\|AP_{\delta}x\|^{2} \\
& = \lim_{\delta\downarrow 0}\int \left|\frac{1+\lambda}{1-\lambda}\right|^{2}|G_{\delta}(\lambda)|^{2}d\|F(\lambda)x\|^{2} \\
& = \int \left|\frac{1+\lambda}{1-\lambda}\right|^{2}d\|F(\lambda)x\|^{2} < \infty.
\end{align}
$$
Conversely if the last integral on the right is finite for some $x$, then the following limit exists in $X$:
$$
y = \lim_{\delta\downarrow 0}\int i\frac{1+\lambda}{1-\lambda}G_{\delta}(\lambda)dF(\lambda)x = \lim_{\delta\downarrow 0}AP_{\delta}x.
$$
Then, because $\lim_{\delta}P_{\delta}x=x$ exists, and because $A$ is closed, it follows that $x\in\mathcal{D}(A)$ and $Ax=y$. Finally, one concludes that
$$
x \in \mathcal{D}(A) \iff \int \left|\frac{1+\lambda}{1-\lambda}\right|^{2}d\|F(\lambda)x\|^{2} < \infty.
$$
And, in that case,
$$
Ax = \lim_{\delta\downarrow 0}\int i \frac{1+\lambda}{1-\lambda}G_{\delta}(\lambda)dF(\lambda)x.
$$
Change of variables: The final step is a change of variables. Define a new spectral measure $E$ on $\mathbb{R}$ by $E(S)=F(\{ \frac{t-i}{t+i} : t\in S\})$. It follows that $\frac{t-i}{t+i}=\lambda$ gives $t=i\frac{1+\lambda}{1-\lambda}$. So,
$$
x \in \mathcal{D}(A) \iff \int_{-\infty}^{\infty}t^{2}d\|E(t)x\|^{2} < \infty,
$$
and, for any such $x$, the following exists as an improper integral:
$$
Ax = \int_{-\infty}^{\infty}tdE(t)x,\;\; x \in \mathcal{D}(A).
$$