Maybe I am missing something but
\begin{equation}
\sum_{n \neq 0} e^{-i E_n (T-i\epsilon)} |n \rangle \langle n|0\rangle =\sum_{n \neq 0} e^{-i E_n T} e^{-\epsilon E_n } |n \rangle \langle n|0\rangle
Β \end{equation}
and not a sum of series. This really allows you then to kill non zero $n$ as $T$ goes to infinity.
As for the physical meaning of this trick, you can imagine it (I imaginie it like that at least) as an account for vanishing mixing instead of purely zero mixing. In other word you add some kind of temperature to your system that you know tend to select out the ground state and then you take the limit when the temperature goes to zero.
Essentially what happens is that the initial state $|0\rangle$ is not an eigenstate of the hamiltonian $H$. Because of that the state evolution is best expressed by expanding it in the eigenstate basis set $|n\rangle$. Now, you can imagine all these states being visited every now and then during the evolution (that's what the above equation states).
A system is always coupled in some way to a thermostat even if its temperature is vanishing. Over a long period, the frequency with which the quantum oscillations will allow each state $|n\rangle$ to be visited will be roughly given by a Boltzmann weight. That's what I meant with the term "mixing", your initial state is "shared" among all the eigenstates and eventually, one can estimate the fraction of it that goes to each $|n\rangle$. In the case of a vanishing temperature, the ground state is greatly favored which explains the result.
Instead of a trick, it is actually more a reminder of what should be the true calculation if it were to be done with the density matrix instead of pure states.
Review of the picture transform
So the "picture transform" is a quantum coordinate transform which basically starts from the usual facts of quantum mechanics. These are roughly: systems are described by vectors in a Hilbert space; everything you can observe is modeled by a Hermitian operator on that Hilbert space, with eigenvalues being the actual values you observe; the only prediction of QM in any particular case is an average value of some observable $A$ in some state $\psi$ given by $\langle A \rangle_\psi = \langle \psi|\hat A |\psi\rangle$; some observable corresponding to total energy $\hat H$ also governs the time evolution of the state vectors in the Hilbert space as $i\hbar \partial_t|\psi\rangle = \hat H |\psi\rangle$.
The picture transform says: hey, that above equation is solved by $|\psi\rangle = U(t) |\psi_0\rangle$ for some unitary $U$ given by $\hat H$ as $$i \hbar \partial_t U = \hat H ~U.$$But suppose the Hamiltonian were instead $\hat h$ for some arbitrary Hermitian $\hat h.$ Then we would find $|\psi\rangle = u(t) |\psi_0\rangle$ for some other unitary $u$. Now insert $u u^\dagger = 1$ into the above to find, $$\langle A\rangle_\psi = \langle \psi_0| U^\dagger~u u^\dagger ~\hat A ~u u^\dagger ~U |\psi_0\rangle.$$Defining $A' = u^\dagger \hat A u$ and $|\psi'\rangle = u^\dagger U|\psi_0\rangle$ must therefore preserve all predictions of quantum mechanics, with the only cost being that now we must consider $$i\hbar \partial_t A' = u^\dagger (\hat A \hat h - \hat h \hat A) u=A'h' - h' A',$$ while for our Hamiltonian we find,$$
i\hbar \partial_t|\psi'\rangle = u^\dagger (- \hat h + \hat H) U |\psi_0\rangle = (H' - h')|\psi'\rangle.$$(Both formulas straightforwardly map to the right-hand-sides after strategically inserting $uu^\dagger = 1$ terms.)
Just like you can mentally understand $U$ as being some sort of time-ordered continuous product and write down a helpful mnemonic,$$U(t) = \mathcal T \exp\left(\frac1{i\hbar} \int_{t_0}^t d\tau ~ \hat H(\tau)\right),$$ if you define $\eta' = H' - h'$ you will find a helpful mnemonic that the new unitary evolution operator looks like,$$U'(t) = u^\dagger U = \mathcal T \exp\left(\frac1{i\hbar} \int_{t_0}^t d\tau ~ \eta'(\tau)\right).$$I claim that this last expression is exactly what was intended in the expression you have above, with $H_I(t)$ taking the place of $\eta'.$
Why the S-matrix tells us everything
The S-matrix is just a matrix expansion of $U'.$ (Equivalently, $U$, when taking the picture transform as the identity transform, the so-called "SchrΓΆdinger picture.") But suppose we have a basis $|n\rangle$ for the Hilbert space, $1 = \sum_n |n\rangle\langle n|.$ Then our above expression for the expectation value can be calculated as:$$\langle A\rangle_\psi = \sum_{mnpq} \langle \psi_0|m\rangle\langle m|{U'}^\dagger|n\rangle\langle n|A'|p\rangle\langle p|U'|q\rangle\langle q|\psi_0\rangle. = (\psi_0)_m^*~U^\dagger_{mn}~A_{np} ~U_{pq}~(\psi_0)_q.$$Therefore $U_{pq} = \langle p | U' | q\rangle$ really does in some sense contain all of the mystery of the time-evolution. And in your case it seems like you're probably dealing with, say, plane waves coming in from infinity getting scattered out to other plane waves travelling out to infinity: this is a great set of states to get a very simple picture, "here is what my Hamiltonian does."
Mixing these together in the interaction picture.
Usually when we're doing this we have some Hamiltonian which we can solve analytically, in your case it seems to be something like $\hat h = \hbar \omega_p a^\dagger_p a_p,$ where $p$ indexes a particle in a definite momentum-state.
When you do this, the presciption above basically just puts phase factors $e^{\pm~i\omega_p t}$ on all of your creation/annihilation operators, and on your definite-momentum states.
If that's right then your question becomes a non-issue: if you paid very careful attention to the $\propto$ sign, and tried to convert it into an $=$ sign, you would have to insert phase factors all over the place in order to keep this thing sane. I think Peskin and Schroeder are probably implicitly saying "hey, you've got this big complicated phase factor out front which does technically matter, but we are explicitly not considering it as important for the physics that we are about to discuss."
I don't have the book but that's where my gut would point me. They are proportional except for some sort of awkward phase factor which would make the equation too long to fit cleanly on one line, so it got absorbed into a $\propto$ symbol.
Best Answer
This is correct. Usually you see it used the opposite way, i.e. $x + X \approx X$ when $X>>x$ but here we are just doing $T \approx T +t_0$ since $T>>T_0$.
No, you can add this because $H_0|0\rangle = 0$, i.e. $e^{xH_0}|0\rangle = |0\rangle .$
We can't use $U(t, t_0) = e^{iH_0(t-t_0)}e^{-iH(t-t_0)}$ here, it simply doesn't apply because we can't change $t_0$. Use Eq. (4.25) where we can set both times freely
$$U(t,t') = e^{iH_0(t-t_0)}e^{-iH(t-t')}e^{-iH_0(t'-t_0)}$$ so for $$U(t_0, -T) = e^{iH_0(t_0-t_0)}e^{-iH(t_0-(-T))}e^{-iH_0((-T)-t_0)}$$ and we obtain whats needed.