Linear response correlation functions I'm working my way through Methods of Molecular Quantum Mechanics by R. McWeeny and have run into a derivation I can't seem to figure out. 
So in chapter 12, he obtains an expression for the first order coefficients of the perturbed wavefunction with respect to a perturbation $H'(t)=F(t)\mathbf{A}$. $\mathbf{A}$ is a hermitian operator and $F(t)$ is a time dependent strength factor and the system was assumed to have started in the state $|0\rangle$ with the perturbation being weak so that these coefficients vary slowly.
$$c_n^{(1)}=(i\hbar)^{-1}\int_{-\infty}^t\langle n|\mathbf{A}|0 \rangle F(t')\exp(i\omega_{n0}t')dt'$$
I'm fine with this expression. Where I get confused is when we try to use this expression to determine the response of some operator $\mathbf{B}$ to the perturbation described by $\mathbf{A}$. He writes that 
$$\langle \mathbf{B} \rangle-\langle \mathbf{B} \rangle_0=\delta\langle \mathbf{B} \rangle=$$
$$(i\hbar)^{-1}\int_{-\infty}^t\sum_{n\neq0}\bigr[ \langle 0|\mathbf{B}|n \rangle \langle n|\mathbf{A}|0 \rangle \exp(-i\omega_{n0}(t-t'))-\langle 0|\mathbf{A}|n \rangle \langle n|\mathbf{B}|0 \rangle \exp(i\omega_{n0}(t-t'))\bigr]F(t')dt'$$
I can't seem to figure out he gets this expression. My thought is to expand $$\langle \Psi'|\delta\mathbf{B}|\Psi' \rangle$$
where $$|\Psi' \rangle=\sum_{n=0} c_n(t)e^{-i\omega_{n0}t}|n\rangle$$
I would hope this would lead to terms like $\langle0|\mathbf{B}|\Psi'\rangle$, but I'm getting extra terms that I can't figure out how to remove.
 A: Try rewriting $\delta\langle \mathbf{B} \rangle= \langle \mathbf{B} \rangle-\langle \mathbf{B} \rangle_0$ more explicitly.
$$\begin{align} \langle \mathbf{B} \rangle-\langle \mathbf{B} \rangle_0 &= \text{Tr}[\mathbf{P}(t)\cdot \mathbf{B}] - \text{Tr}[\mathbf{P}(0)\cdot \mathbf{B}] \\ &= \langle 0 | \mathbf{P}(t) \mathbf{B} | 0 \rangle - \langle 0 | \mathbf{P}(0) \mathbf{B} | 0 \rangle  \\  \end{align}$$
Recall that to first order, the density matrix $\mathbf{P}(t)$ can be written
$$\mathbf{P}(t) = |0\rangle\langle0| + |\Psi^{(1)}(t)\rangle\langle0| + |0\rangle \langle \Psi^{(1)}(t)| + \cdots $$
You can plug in your expression
$$|\Psi' \rangle = |\Psi^{(1)}(t)\rangle =\sum_{n=0} c^{(1)}_n(t)e^{-i\omega_{n0}t}|n\rangle$$
This should give you (treating $\mathbf{B}$ as Hermitian operator)
$$\delta\langle \mathbf{B} \rangle = \langle 0 | \mathbf{B}\sum_{n=0} c^{(1)}_n(t)e^{-i\omega_{n0}t}|n\rangle\langle0|0\rangle + \langle 0 | 0 \rangle  \sum_{n=0} 
 \langle n |  c^{(1)*}_n(t)e^{+i\omega_{n0}t}\mathbf{B} |0\rangle $$
From here it is a matter of plugging in your first expression
$$c_n^{(1)}=(i\hbar)^{-1}\int_{-\infty}^t\langle n|\mathbf{A}|0 \rangle F(t')\exp(i\omega_{n0}t')dt'$$
And cleaning up until you get 
$$(i\hbar)^{-1}\int_{-\infty}^t\sum_{n\neq0}\bigr[ \langle 0|\mathbf{B}|n \rangle \langle n|\mathbf{A}|0 \rangle \exp(-i\omega_{n0}(t-t'))-\langle 0|\mathbf{A}|n \rangle \langle n|\mathbf{B}|0 \rangle \exp(i\omega_{n0}(t-t'))\bigr]F(t')dt'$$
