Recovering phases of a sum of sinusoids using a DFT

I generate data in the following way: \begin{align} x_n = \cos\left(\left(\frac{2\pi}{100}\cdot 1\right)n + 0.8\right) + \cos\left(\left(\frac{2\pi}{100}\cdot 2\right)n + 0.6\right) + \cos\left(\left(\frac{2\pi}{100}\cdot 3\right)n + 0.4\right) + \cos\left(\left(\frac{2\pi}{100}\cdot 4\right)n + 0.2\right) \end{align} I want to use the DFT to recover the phases 0.8,0.6,0.4,0.2. Here is my code (written in R):

q = 100
ts = 1:q
x = cos((2*pi*1/q)*ts + 1) + cos((2*pi*2/q)*ts + 0.8) +
cos((2*pi*3/q)*ts + 0.6) + cos((2*pi*4/q)*ts + 0.4)
plot(ts, x)
X <- fft(x)


I would have thought the way to recover the phases would be to just use $$\phi_k = \text{atan2}(\text{Im}(X_k),\text{Re}(X_k))$$, but it doesn't appear to give the correct values.

phases = rep(0, 4)
for(k in 1:4) {
phases[k] <- atan2(Im(X[k+1]), Re(X[k+1]))
}
print(phases)
[1] 1.0628319 0.9256637 0.7884956 0.6513274


However, if I let $$\phi_k' = \phi_k - \frac{2\pi k}{100}$$, then I get the correct values. Why is this the case?

for(k in 1:4) {
phases[k] <- phases[k] - 2*pi*k/q
}
print(phases)
[1] 1.0 0.8 0.6 0.4


Because you start your time vector at 1 and not at zero. You should use ts = 0:(q-1). Starting at one is equivalent to a delay of -1 sample delay the transfer function of which is $$e^{j 2\pi \omega*T}$$ where T is your sample period. For the FFT that comes out to be $$e^{j 2\pi k/100}$$