Entirely linear neural network learning non-linear function

Question

I have a neural network that's trained on a sine wave. It uses a lookback of 20 to see what the last 20 predictions were and predict the next value. This network has only a single Linear layer (input size 20, output size 1) with no activation function and from just that, it is able to extrapolate a sine wave almost perfectly.

This is very confusing to me, as the $sin$ function is non-linear: $sin(a+b) \neq sin(a) + sin(b)$ so the network shouldn't be able to approximate it (very well).

The code to reproduce this is below:

import numpy as np
import matplotlib.pyplot as plt

import torch
import torch.nn as nn
from torch.utils.data import DataLoader, Dataset
import torch.optim as optim

X_train = np.arange(0,100,0.5) 
y_train = np.sin(X_train)

X_test = np.arange(100,200,0.5) 
y_test = np.sin(X_test)

n_features = 1

train_series = torch.from_numpy(y_train)
test_series = torch.from_numpy(y_test)

# Expects input of (batch, sequence, features)
# So shape should be (1, 179, 20) and labels (1, 1, 179)
look_back = 20

train_dataset = []
train_labels = []
for i in range(len(train_series)-look_back):
    train_dataset.append(train_series[i:i+20])
    train_labels.append(train_series[i+20])
train_dataset = torch.stack(train_dataset).unsqueeze(0)
train_labels = torch.stack(train_labels).unsqueeze(0).unsqueeze(2)

class Net(nn.Module):
    def __init__(self, input_shape):
        super(Net, self).__init__()
        
        self.fc = nn.Linear(input_shape, 1)
    
    def forward(self, x):
        out = self.fc(x)
        return out


model = Net(look_back).double()
loss_function = nn.MSELoss()
optimizer = optim.Adam(model.parameters(), lr=0.001)

loss_curve = []
for epoch in range(300):
    loss_total = 0
    
    model.zero_grad()
    
    predictions = model(train_dataset)
    
    loss = loss_function(predictions, train_labels)
    loss_total += loss.item()
    loss.backward()
    optimizer.step()
    loss_curve.append(loss_total)

extrapolation = []
seed_batch = test_series[:20].reshape(1, 1, 20)
current_batch = seed_batch
with torch.no_grad():
    for i in range(180):
        predicted_value = model(current_batch)
        extrapolation.append(predicted_value.item())
        current_batch = torch.cat((current_batch[:,:,1:], predicted_value), axis=2)

x = np.arange(110,200,0.5)
fig, ax = plt.subplots(1, 1, figsize=(15, 5))
ax.plot(X_train,y_train, lw=2, label='train data')
ax.plot(X_test,y_test, lw=3, c='y', label='test data')
ax.plot(x,extrapolation, lw=3, c='r',linestyle = ':', label='extrapolation')
ax.legend(loc="lower left")
plt.show();

And it produces the following plot:

Answer 1

You can even predict a sinusoid with much less than 20 samples: two previous samples suffice. The reason is that sinusoids appear as solutions of second-order linear difference equations.

In other words, you can write the linear relationship of your sinusoid with frequency $f$ and phase $\phi$ as

$sin (f \times n + \phi) = a_1 sin (f \times (n-1) + \phi) + a_2 sin (f \times (n-2) + \phi)$

with real prediction coefficients $a_1$ and $a_2$ and discrete time index $n$.