<jupyter_start><jupyter_text>Distributed lag model$$y_{t}=a+w_{0}x_{t}+w_{1}x_{t-1}+w_{2}x_{t-2}+...+{\text{error term}}$$https://en.wikipedia.org/wiki/Distributed_lag<jupyter_code>import numpy as np
import matplotlib.pyplot as plt

N = 100
time = np.arange(N+1)
dt = 1
d = 3

x_diff = np.random.normal(size = [N,d]) * np.sqrt(dt)
x = np.cumsum(x_diff,axis = 0)
x = np.concatenate([np.zeros([1,d]), x], axis = 0)

# x1 = 0.5 * np.sqrt(time)[:,None]
# x2 = np.sin(time)[:,None]
# x3 = np.log(time+1)[:,None]
# x = np.concatenate([x1,x2,x3], axis = 1)

plt.plot(time,x)
plt.show()
def lag(x,length,gap):
    x = x[:,None,:]
    X = x[:1-length-gap]
    for i in range(1,length):
        helper = x[i:1-length-gap+i]
        X = np.concatenate([X, helper], axis = 1)
    return X
x_lag = lag(x,3,1)
len_lag = x_lag.shape[0]
import torch
import torch.nn as nn
x_lag_flatten = torch.tensor(x_lag, dtype = torch.float).view([len_lag,-1])
linear = nn.Linear(x_lag_flatten.shape[-1],1)
y = linear(x_lag_flatten)
y = y.detach().numpy() + np.random.normal(size = [len_lag,1])
plt.plot(y)
plt.show()
import sklearn
split = 60
num_sample = x_lag.shape[0]
X = np.reshape(x_lag,[num_sample,-1])
y = y
X_train = X[:split]
y_train = y[:split]
X_test = X[split:]
y_test = y[split:]
from sklearn.linear_model import Ridge
model = Ridge(alpha=1)
model.fit(X_train, y_train)
y_train_predict = model.predict(X_train)
y_test_predict = model.predict(X_test)
y_predict = np.concatenate([y_train_predict,y_test_predict],axis = 0)
plt.plot(y)
plt.plot(y_predict)
lmax = max(y.max(),y_predict.max())
lmin = min(y.min(),y_predict.min())
plt.vlines(split,lmin,lmax,'r')
plt.legend(['True', 'Predict', 'train-test split'])
plt.show()
def lag(x,length,gap):
    x = x[:,None,:]
    X = x[:1-length-gap]
    for i in range(1,length):
        helper = x[i:1-length-gap+i]
        X = np.concatenate([X, helper], axis = 1)
    return X<jupyter_output><empty_output><jupyter_text>Real Data $$d_{t} = f(c_{t-1}, c_{t-2}, c_{t-3}) + noise$$where $d$ is death and $c$ is new cases<jupyter_code>import pandas as pd
from sklearn import preprocessing
import matplotlib.pyplot as plt
from sklearn.linear_model import Ridge
df = pd.read_csv("data/1_owid/owid-covid-data.csv")
df_europe = df[df['continent'] == 'Europe']
df0 = df_europe
df0 = df0[['date',\
           'location',\
       'new_cases_smoothed_per_million',\
       'new_deaths_smoothed_per_million',\
       'icu_patients_per_million']]
#        'positive_rate']]
df0 =df0.dropna()
df00 = df0

C = list(df0['location'].drop_duplicates())
df00 = df0.drop(columns=['date','location'])
x = df00.values 
min_max_scaler = preprocessing.MinMaxScaler()
x_scaled = min_max_scaler.fit_transform(x)
df00 = pd.DataFrame(x_scaled, index = df00.index, columns=df00.columns)

df_helper = df_europe[['location', 'population','population_density']]
POPULATION = []
POPULATION_DENSITY = []
for c in C:
    population = df_helper[df_helper['location'] == c].iloc[0]['population']
    density = df_helper[df_helper['location'] == c].iloc[0]['population_density']
    POPULATION.append(population)
    POPULATION_DENSITY.append(density)
C_large = []
POPULATION_large = []
choose_large = []
for c,p in zip(C,POPULATION):
    if p > 9*10**6:
        C_large.append(c)
        POPULATION_large.append(p)
C_large<jupyter_output><empty_output><jupyter_text>Prediction<jupyter_code>DICT = {}
TRAIN = []
TEST = []
RANGE_TRAIN = []
RANGE_TEST = []
for c in C_large:
    print(c)
    cc = df00[df0['location'] == c]
    x = cc[['new_cases_smoothed_per_million', 'icu_patients_per_million', 'new_deaths_smoothed_per_million']].values
    y = cc['new_deaths_smoothed_per_million'].values
    length = 5
    gap = 5
    x_lag = lag(x,length,gap)
    y_lag = y[:len(x_lag),None]
    
    num_sample = x_lag.shape[0]
    split = int(0.7*num_sample)
    X = np.reshape(x_lag,[num_sample,-1])
    Y = y_lag
    X_train = X[:split]
    Y_train = Y[:split]
    X_test = X[split:]
    Y_test = Y[split:]
    
    model = Ridge(alpha=1)
    model.fit(X_train, Y_train)
    
    train_score = model.score(X_train,Y_train)
    test_score = model.score(X_test,Y_test)
    print('train', train_score)
    print('test', test_score)
    TRAIN.append(train_score)
    TEST.append(test_score)
    
    DICT[c] = [cc, model, X_train, Y_train, X_test, Y_test, X, Y]
    
    Y_train_predict = model.predict(X_train)
    Y_test_predict = model.predict(X_test)
    Y_predict = np.concatenate([Y_train_predict,Y_test_predict],axis = 0)
    Y_predict_aug = np.concatenate([Y[:length+gap-1,:],Y_predict],axis = 0)
    
    cc = cc.assign(predict_death=pd.Series(Y_predict_aug[:,0]).values)
    
    range_train = [x_lag[:split,:,1].min(),x_lag[:split,:,1].max()]
    range_test = [x_lag[split:,:,1].min(),x_lag[split:,:,1].max()]
    print('range_train: ', range_train)
    print('range_test : ', range_test )
    RANGE_TRAIN.append(range_train )
    RANGE_TEST.append(range_test )
    cc.plot()
    
#     plt.plot(Y)
#     plt.plot(Y_predict)
#     lmax = max(y.max(),Y_predict.max())
#     lmin = min(y.min(),Y_predict.min())
    plt.vlines(cc.index[split],-0.1,1.1,'r')
#     plt.legend(['True', 'Predict', 'train-test split'])
    plt.show()
    
import string
from matplotlib import cm

def plot_scatter_helper(name_x, name_y, x, y, C):
    df_helper = pd.DataFrame({name_x:x, name_y:y}, 
                      index=C)
    fig, ax = plt.subplots(figsize = [12,4])


    cmap = cm.get_cmap('Spectral')
    df_helper.plot(name_x, name_y, kind='scatter', ax=ax, s=100, linewidth=0, 
            c=range(len(df_helper)), colormap=cmap)

    for k, v in df_helper.iterrows():
        ax.annotate(k, v,
                    xytext=(10,-5), textcoords='offset points',
                    family='sans-serif', fontsize=10, color='darkslategrey')
    fig.canvas.draw()<jupyter_output><empty_output><jupyter_text>No relation between population and score<jupyter_code>name_x = 'LogPopulation'
name_y = 'transformed test score'
x = np.log(np.array(POPULATION_large))
y = np.exp(np.array(TEST)-1)
print('LogPopulation - transformed train score')
plot_scatter_helper(name_x, name_y, x, y, C_large)

# name_x = 'LogPopulation'
# name_y = 'transformed train score'
# x = np.log(np.array(POPULATION))
# y = np.exp(np.array(TRAIN)-1)
# print('train score')
# plot_scatter_helper(name_x, name_y, x, y)<jupyter_output>LogPopulation - transformed train score<jupyter_text>Strong Relation between peak and score<jupyter_code>def overlap(train, test):
    error1 = train[0] - test[0]
    error1 = (error1 > 0) * error1
    error2 = test[1] - train[1]
    error2 = (error2 > 0) * error2 
    return error1 + error2
DIFF_lap = [overlap(train,test) for train,test in zip(RANGE_TRAIN,RANGE_TEST)]

DIFF = [train[1] - train[0] for train in RANGE_TRAIN]

name_x = 'Suffering from first round'
name_y = 'transformed test score'
x = np.array(DIFF) - np.array(DIFF_lap)
y = np.exp((np.array(TEST)-1))
print('Suffering from first round  - transformed train score')
plot_scatter_helper(name_x, name_y, x, y, C_large)

model = Ridge(alpha=1)
model.fit(x[:,None], y[:,None])
lin_x = np.linspace(x.min(),x.max(),10)
lin_y = model.predict(lin_x[:,None])
plt.plot(lin_x,lin_y,alpha = 0.5)
# plt.plot([x.min(),x.max()],[y.min(),y.max()],alpha = 0.5)
plt.legend(['Linear regression'])
plt.show()<jupyter_output>Suffering from first round  - transformed train score<jupyter_text>Use Italy to predict Austria<jupyter_code>_, model, _, _, _, _, _, _ =  DICT['Italy']
cc, model0, X_train, Y_train, X_test, Y_test, X, Y = DICT['Austria']

Y_train_predict = model.predict(X_train)
Y_test_predict = model.predict(X_test)
Y_predict = np.concatenate([Y_train_predict,Y_test_predict],axis = 0)
Y_predict_aug = np.concatenate([Y[:length+gap-1,:],Y_predict],axis = 0)

Y_train_predict0 = model0.predict(X_train)
Y_test_predict0 = model0.predict(X_test)
Y_predict0 = np.concatenate([Y_train_predict0,Y_test_predict0],axis = 0)
Y_predict_aug0 = np.concatenate([Y[:length+gap-1,:],Y_predict0],axis = 0)

cc = cc.assign(predict_death_Italy=pd.Series(Y_predict_aug[:,0]).values)
cc = cc.assign(predict_death_Austria=pd.Series(Y_predict_aug0[:,0]).values)
cc.plot(figsize = [12,4])
#     plt.plot(Y)
#     plt.plot(Y_predict)
#     lmax = max(y.max(),Y_predict.max())
#     lmin = min(y.min(),Y_predict.min())

num_samples = X.shape[0]
split = int(0.7*num_samples )
plt.vlines(cc.index[split],-0.1,1.1,'r')
#     plt.legend(['True', 'Predict', 'train-test split'])
plt.show()<jupyter_output><empty_output><jupyter_text>Positive rate NOT helps<jupyter_code>df = pd.read_csv("data/1_owid/owid-covid-data.csv")
df_europe = df[df['continent'] == 'Europe']
df0 = df_europe
df0 = df0[['date',\
           'location',\
       'new_cases_smoothed_per_million',\
       'new_deaths_smoothed_per_million',\
       'icu_patients_per_million',\
       'positive_rate']]
df0 =df0.dropna()
df00 = df0

C = list(df0['location'].drop_duplicates())
df00 = df0.drop(columns=['date','location'])
x = df00.values 
min_max_scaler = preprocessing.MinMaxScaler()
x_scaled = min_max_scaler.fit_transform(x)
df00 = pd.DataFrame(x_scaled, index = df00.index, columns=df00.columns)

df_helper = df_europe[['location', 'population','population_density']]
POPULATION = []
POPULATION_DENSITY = []
for c in C:
    population = df_helper[df_helper['location'] == c].iloc[0]['population']
    density = df_helper[df_helper['location'] == c].iloc[0]['population_density']
    POPULATION.append(population)
    POPULATION_DENSITY.append(density)
    
C_large = []
POPULATION_large = []
choose_large = []
for c,p in zip(C,POPULATION):
    if p > 9*10**6:
        C_large.append(c)
        POPULATION_large.append(p)
C_large

DICT = {}
TRAIN = []
TEST_positive = []
RANGE_TRAIN = []
RANGE_TEST = []
for c in C_large:
#     print(c)
    cc = df00[df0['location'] == c]
    x = cc[['new_cases_smoothed_per_million', 'icu_patients_per_million', 'new_deaths_smoothed_per_million', 'positive_rate']].values
    y = cc['new_deaths_smoothed_per_million'].values
    length = 5
    gap = 5
    x_lag = lag(x,length,gap)
    y_lag = y[:len(x_lag),None]
    
    num_sample = x_lag.shape[0]
    split = int(0.7*num_sample)
    X = np.reshape(x_lag,[num_sample,-1])
    Y = y_lag
    X_train = X[:split]
    Y_train = Y[:split]
    X_test = X[split:]
    Y_test = Y[split:]
    
    model = Ridge(alpha=1)
    model.fit(X_train, Y_train)
    
    train_score = model.score(X_train,Y_train)
    test_score = model.score(X_test,Y_test)
#     print('train', train_score)
#     print('test', test_score)
    TRAIN.append(train_score)
    TEST_positive.append(test_score)
    
    DICT[c] = [cc, model, X_train, Y_train, X_test, Y_test, X, Y]
    
    Y_train_predict = model.predict(X_train)
    Y_test_predict = model.predict(X_test)
    Y_predict = np.concatenate([Y_train_predict,Y_test_predict],axis = 0)
    Y_predict_aug = np.concatenate([Y[:length+gap-1,:],Y_predict],axis = 0)
    
    cc = cc.assign(predict_death=pd.Series(Y_predict_aug[:,0]).values)
    
    range_train = [x_lag[:split,:,1].min(),x_lag[:split,:,1].max()]
    range_test = [x_lag[split:,:,1].min(),x_lag[split:,:,1].max()]
#     print('range_train: ', range_train)
#     print('range_test : ', range_test )
    RANGE_TRAIN.append(range_train )
    RANGE_TEST.append(range_test )
#     cc.plot()
    
# #     plt.plot(Y)
# #     plt.plot(Y_predict)
# #     lmax = max(y.max(),Y_predict.max())
# #     lmin = min(y.min(),Y_predict.min())
#     plt.vlines(cc.index[split],-0.1,1.1,'r')
# #     plt.legend(['True', 'Predict', 'train-test split'])
#     plt.show()
    
name_x = 'TEST'
name_y = 'TEST_use_positive_rate'
x = np.exp(np.array(TEST)-1)
y = np.exp(np.array(TEST_positive)-1)
print('TEST - TEST_use_positive_rate')
plot_scatter_helper(name_x, name_y, x, y, C_large)
plt.plot([0,1],[0,1],alpha = 0.5)
plt.legend(['y = x'])
plt.show()<jupyter_output>TEST - TEST_use_positive_rate<jupyter_text>Therefore accuracy no longer depends on existence of peak in cases because positive rate captures deaths<jupyter_code>name_x = 'Suffering from first round'
name_y = 'transformed test score with positive'
x = np.array(DIFF) - np.array(DIFF_lap)
y = np.exp((np.array(TEST_positive)-1))
print('Suffering from first round  - transformed test score with positive')
plot_scatter_helper(name_x, name_y, x, y, C_large)<jupyter_output>Suffering from first round  - transformed test score with positive
{'fraction_non_alphanumeric': 0.16670689516453455, 'fraction_numerical': 0.018022367044814545, 'mean_word_length': 4.444590451160753, 'pattern_counts': {'":': 0, '<': 29, '<?xml version=': 0, '>': 33, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '311706', 'license': 'Apache-2.0', 'n_tokens_mistral': 5161, 'n_tokens_neox': 4984, 'n_words': 810, 'path': 'distributed lag model.ipynb', 'repo_name': 'justinhou95/CovidAnalysis', 'size': '12429'}