text_plain_output_1.png
from sklearn.metrics import mean_squared_error

from sklearn.metrics import mean_squared_error, mean_squared_log_error

import numpy as np

import numpy as np # linear algebra

import pandas as pd

import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)

import xgboost as xgb

class CovidModel:
    def __init__(self):
        pass
    def predict_first_day(self, date):
        return None
    def predict_next_day(self, yesterday_pred_df):
        return None
class CovidModelGIBA(CovidModel):
    def __init__(self, lag=1, seed=1):
        self.lag = lag
        self.seed = seed
        train = pd.read_csv('../input/covid19-global-forecasting-week-4/train.csv')
        train['Date'] = pd.to_datetime(train['Date'])
        self.maxdate = str(train['Date'].max())[:10]
        self.testdate = str(train['Date'].max() + pd.Timedelta(days=1))[:10]
        train['Province_State'].fillna('', inplace=True)
        train['day'] = train.Date.dt.dayofyear
        self.day_min = train['day'].min()
        train['day'] -= self.day_min
        train['geo'] = ['_'.join(x) for x in zip(train['Country_Region'], train['Province_State'])]
        test = pd.read_csv('../input/covid19-global-forecasting-week-4/test.csv')
        test['Date'] = pd.to_datetime(test['Date'])
        test['Province_State'].fillna('', inplace=True)
        test['day'] = test.Date.dt.dayofyear
        test['day'] -= self.day_min
        test['geo'] = ['_'.join(x) for x in zip(test['Country_Region'], test['Province_State'])]
        test['Id'] = -1
        test['ConfirmedCases'] = 0
        test['Fatalities'] = 0
        self.trainmaxday = train['day'].max()
        self.testday1 = train['day'].max() + 1
        self.testdayN = test['day'].max()
        publictest = test.loc[test.Date > train.Date.max()].copy()
        train = pd.concat((train, publictest), sort=False)
        train.sort_values(['Country_Region', 'Province_State', 'Date'], inplace=True)
        train = train.reset_index(drop=True)
        train['ForecastId'] = pd.merge(train, test, on=['Country_Region', 'Province_State', 'Date'], how='left')['ForecastId_y'].values
        train['cid'] = train['Country_Region'] + '_' + train['Province_State']
        train['log0'] = np.log1p(train['ConfirmedCases'])
        train['log1'] = np.log1p(train['Fatalities'])
        train = train.loc[(train.log0 > 0) | train.ForecastId.notnull() | (train.Date >= '2020-03-17')].copy()
        train = train.reset_index(drop=True)
        train['days_since_1case'] = train.groupby('cid')['Id'].cumcount()
        dt = pd.read_csv('../input/covid19-lockdown-dates-by-country/countryLockdowndates.csv')
        dt.columns = ['Country_Region', 'Province_State', 'Date', 'Type', 'Reference']
        dt = dt.loc[dt.Date == dt.Date]
        dt['Province_State'] = dt['Province_State'].fillna('')
        dt['Date'] = pd.to_datetime(dt['Date'])
        dt['Date'] = dt['Date'] + pd.Timedelta(days=8)
        dt['Type'] = pd.factorize(dt['Type'])[0]
        dt['cid'] = dt['Country_Region'] + '_' + dt['Province_State']
        del dt['Reference'], dt['Country_Region'], dt['Province_State']
        train = pd.merge(train, dt, on=['cid', 'Date'], how='left')
        train['Type'] = train.groupby('cid')['Type'].fillna(method='ffill')
        train['target0'] = np.log1p(train['ConfirmedCases'])
        train['target1'] = np.log1p(train['Fatalities'])
        self.train = train.copy()
    def create_features(self, df, valid_day):
        df = df.loc[df.day >= valid_day - 50].copy()
        df['lag0_1'] = df.groupby('cid')['target0'].shift(self.lag)
        df['lag0_1'] = df.groupby('cid')['lag0_1'].fillna(method='bfill')
        df['lag0_8'] = df.groupby('cid')['target0'].shift(8)
        df['lag0_8'] = df.groupby('cid')['lag0_8'].fillna(method='bfill')
        df['lag1_1'] = df.groupby('cid')['target1'].shift(self.lag)
        df['lag1_1'] = df.groupby('cid')['lag1_1'].fillna(method='bfill')
        df['m0'] = df.groupby('cid')['lag0_1'].rolling(2).mean().values
        df['m1'] = df.groupby('cid')['lag0_1'].rolling(3).mean().values
        df['m2'] = df.groupby('cid')['lag0_1'].rolling(4).mean().values
        df['m3'] = df.groupby('cid')['lag0_1'].rolling(5).mean().values
        df['m4'] = df.groupby('cid')['lag0_1'].rolling(7).mean().values
        df['m5'] = df.groupby('cid')['lag0_1'].rolling(10).mean().values
        df['m6'] = df.groupby('cid')['lag0_1'].rolling(12).mean().values
        df['m7'] = df.groupby('cid')['lag0_1'].rolling(16).mean().values
        df['m8'] = df.groupby('cid')['lag0_1'].rolling(20).mean().values
        df['m9'] = df.groupby('cid')['lag0_1'].rolling(25).mean().values
        df['n0'] = df.groupby('cid')['lag1_1'].rolling(2).mean().values
        df['n1'] = df.groupby('cid')['lag1_1'].rolling(3).mean().values
        df['n2'] = df.groupby('cid')['lag1_1'].rolling(4).mean().values
        df['n3'] = df.groupby('cid')['lag1_1'].rolling(5).mean().values
        df['n4'] = df.groupby('cid')['lag1_1'].rolling(7).mean().values
        df['n5'] = df.groupby('cid')['lag1_1'].rolling(10).mean().values
        df['n6'] = df.groupby('cid')['lag1_1'].rolling(12).mean().values
        df['n7'] = df.groupby('cid')['lag1_1'].rolling(16).mean().values
        df['n8'] = df.groupby('cid')['lag1_1'].rolling(20).mean().values
        df['m0'] = df.groupby('cid')['m0'].fillna(method='bfill')
        df['m1'] = df.groupby('cid')['m1'].fillna(method='bfill')
        df['m2'] = df.groupby('cid')['m2'].fillna(method='bfill')
        df['m3'] = df.groupby('cid')['m3'].fillna(method='bfill')
        df['m4'] = df.groupby('cid')['m4'].fillna(method='bfill')
        df['m5'] = df.groupby('cid')['m5'].fillna(method='bfill')
        df['m6'] = df.groupby('cid')['m6'].fillna(method='bfill')
        df['m7'] = df.groupby('cid')['m7'].fillna(method='bfill')
        df['m8'] = df.groupby('cid')['m8'].fillna(method='bfill')
        df['m9'] = df.groupby('cid')['m9'].fillna(method='bfill')
        df['n0'] = df.groupby('cid')['n0'].fillna(method='bfill')
        df['n1'] = df.groupby('cid')['n1'].fillna(method='bfill')
        df['n2'] = df.groupby('cid')['n2'].fillna(method='bfill')
        df['n3'] = df.groupby('cid')['n3'].fillna(method='bfill')
        df['n4'] = df.groupby('cid')['n4'].fillna(method='bfill')
        df['n5'] = df.groupby('cid')['n5'].fillna(method='bfill')
        df['n6'] = df.groupby('cid')['n6'].fillna(method='bfill')
        df['n7'] = df.groupby('cid')['n7'].fillna(method='bfill')
        df['n8'] = df.groupby('cid')['n8'].fillna(method='bfill')
        df['flag_China'] = 1 * (df['Country_Region'] == 'China')
        df['flag_US'] = 1 * (df['Country_Region'] == 'US')
        df['flag_Kosovo_'] = 1 * (df['cid'] == 'Kosovo_')
        df['flag_Korea'] = 1 * (df['cid'] == 'Korea, South_')
        df['flag_Nepal_'] = 1 * (df['cid'] == 'Nepal_')
        df['flag_Holy See_'] = 1 * (df['cid'] == 'Holy See_')
        df['flag_Suriname_'] = 1 * (df['cid'] == 'Suriname_')
        df['flag_Ghana_'] = 1 * (df['cid'] == 'Ghana_')
        df['flag_Togo_'] = 1 * (df['cid'] == 'Togo_')
        df['flag_Malaysia_'] = 1 * (df['cid'] == 'Malaysia_')
        df['flag_US_Rhode'] = 1 * (df['cid'] == 'US_Rhode Island')
        df['flag_Bolivia_'] = 1 * (df['cid'] == 'Bolivia_')
        df['flag_China_Tib'] = 1 * (df['cid'] == 'China_Tibet')
        df['flag_Bahrain_'] = 1 * (df['cid'] == 'Bahrain_')
        df['flag_Honduras_'] = 1 * (df['cid'] == 'Honduras_')
        df['flag_Bangladesh'] = 1 * (df['cid'] == 'Bangladesh_')
        df['flag_Paraguay_'] = 1 * (df['cid'] == 'Paraguay_')
        tr = df.loc[df.day < valid_day].copy()
        vl = df.loc[df.day == valid_day].copy()
        tr = tr.loc[tr.lag0_1 > 0].copy()
        maptarget0 = tr.groupby('cid')['target0'].agg(log0_max='max').reset_index()
        maptarget1 = tr.groupby('cid')['target1'].agg(log1_max='max').reset_index()
        vl['log0_max'] = pd.merge(vl, maptarget0, on='cid', how='left')['log0_max'].values
        vl['log1_max'] = pd.merge(vl, maptarget1, on='cid', how='left')['log1_max'].values
        vl['log0_max'] = vl['log0_max'].fillna(0)
        vl['log1_max'] = vl['log1_max'].fillna(0)
        return (tr, vl)
    def train_models(self, valid_day=10):
        train = self.train.copy()
        train.loc[(train.cid == 'China_Guizhou') & (train.Date == '2020-03-17'), 'target0'] = np.log1p(146)
        train.loc[(train.cid == 'Guyana_') & (train.Date >= '2020-03-22') & (train.Date <= '2020-03-30'), 'target0'] = np.log1p(12)
        train.loc[(train.cid == 'US_Virgin Islands') & (train.Date >= '2020-03-29') & (train.Date <= '2020-03-29'), 'target0'] = np.log1p(24)
        train.loc[(train.cid == 'US_Virgin Islands') & (train.Date >= '2020-03-30') & (train.Date <= '2020-03-30'), 'target0'] = np.log1p(27)
        train.loc[(train.cid == 'Iceland_') & (train.Date >= '2020-03-15') & (train.Date <= '2020-03-15'), 'target1'] = np.log1p(0)
        train.loc[(train.cid == 'Kazakhstan_') & (train.Date >= '2020-03-20') & (train.Date <= '2020-03-20'), 'target1'] = np.log1p(0)
        train.loc[(train.cid == 'Serbia_') & (train.Date >= '2020-03-26') & (train.Date <= '2020-03-26'), 'target1'] = np.log1p(5)
        train.loc[(train.cid == 'Serbia_') & (train.Date >= '2020-03-27') & (train.Date <= '2020-03-27'), 'target1'] = np.log1p(6)
        train.loc[(train.cid == 'Slovakia_') & (train.Date >= '2020-03-22') & (train.Date <= '2020-03-31'), 'target1'] = np.log1p(1)
        train.loc[(train.cid == 'US_Hawaii') & (train.Date >= '2020-03-25') & (train.Date <= '2020-03-31'), 'target1'] = np.log1p(1)
        param = {'subsample': 1.0, 'colsample_bytree': 0.85, 'max_depth': 5, 'gamma': 0.0, 'learning_rate': 0.01, 'min_child_weight': 6.0, 'reg_alpha': 0.0, 'reg_lambda': 0.4, 'silent': 1, 'objective': 'reg:squarederror', 'nthread': 12, 'seed': self.seed}
        tr, vl = self.create_features(train.copy(), valid_day)
        features = [f for f in tr.columns if f not in ['lag0_8', 'Id', 'ConfirmedCases', 'Fatalities', 'log0', 'log1', 'target0', 'target1', 'ypred0', 'ypred1', 'Province_State', 'Country_Region', 'Date', 'ForecastId', 'cid', 'geo', 'day', 'GDP_region', 'TRUE POPULATION', 'pct_in_largest_city', ' TFR ', ' Avg_age ', 'latitude', 'longitude', 'abs_latitude', 'temperature', 'humidity', 'Personality_pdi', 'Personality_idv', 'Personality_mas', 'Personality_uai', 'Personality_ltowvs', 'Personality_assertive', 'personality_perform', 'personality_agreeableness', 'murder', 'High_rises', 'max_high_rises', 'AIR_CITIES', 'AIR_AVG', 'continent_gdp_pc', 'continent_happiness', 'continent_generosity', 'continent_corruption', 'continent_Life_expectancy']]
        self.features0 = features
        features = [f for f in tr.columns if f not in ['m0', 'm1', 'm2', 'm3', 'Id', 'ConfirmedCases', 'Fatalities', 'log0', 'log1', 'target0', 'target1', 'ypred0', 'ypred1', 'Province_State', 'Country_Region', 'Date', 'ForecastId', 'cid', 'geo', 'day', 'GDP_region', 'TRUE POPULATION', 'pct_in_largest_city', ' TFR ', ' Avg_age ', 'latitude', 'longitude', 'abs_latitude', 'temperature', 'humidity', 'Personality_pdi', 'Personality_idv', 'Personality_mas', 'Personality_uai', 'Personality_ltowvs', 'Personality_assertive', 'personality_perform', 'personality_agreeableness', 'murder', 'High_rises', 'max_high_rises', 'AIR_CITIES', 'AIR_AVG', 'continent_gdp_pc', 'continent_happiness', 'continent_generosity', 'continent_corruption', 'continent_Life_expectancy']]
        self.features1 = features
        nrounds0 = 680
        nrounds1 = 630
        dtrain = xgb.DMatrix(tr[self.features0], tr['target0'])
        param['seed'] = self.seed
        self.model0 = xgb.train(param, dtrain, nrounds0, verbose_eval=0)
        param['seed'] = self.seed + 1
        self.model1 = xgb.train(param, dtrain, nrounds0, verbose_eval=0)
        dtrain = xgb.DMatrix(tr[self.features1], tr['target1'])
        param['seed'] = self.seed
        self.model2 = xgb.train(param, dtrain, nrounds1, verbose_eval=0)
        param['seed'] = self.seed + 1
        self.model3 = xgb.train(param, dtrain, nrounds1, verbose_eval=0)
        self.vl = vl
        return 1
    def predict_first_day(self, day):
        self.day = day
        self.train_models(day)
        dvalid = xgb.DMatrix(self.vl[self.features0])
        ypred0 = (self.model0.predict(dvalid) + self.model1.predict(dvalid)) / 2
        dvalid = xgb.DMatrix(self.vl[self.features1])
        ypred1 = (self.model2.predict(dvalid) + self.model3.predict(dvalid)) / 2
        self.vl['ypred0'] = ypred0
        self.vl['ypred1'] = ypred1
        self.vl.loc[self.vl.ypred0 < self.vl.log0_max, 'ypred0'] = self.vl.loc[self.vl.ypred0 < self.vl.log0_max, 'log0_max']
        self.vl.loc[self.vl.ypred1 < self.vl.log1_max, 'ypred1'] = self.vl.loc[self.vl.ypred1 < self.vl.log1_max, 'log1_max']
        VALID = self.vl[['geo', 'day', 'ypred0', 'ypred1']].copy()
        VALID.columns = ['geo', 'day', 'ConfirmedCases', 'Fatalities']
        return VALID.reset_index(drop=True)
    def predict_next_day(self, yesterday):
        self.day += 1
        feats = ['geo', 'day']
        self.train['ypred0'] = pd.merge(self.train[feats], yesterday[feats + ['ConfirmedCases']], on=feats, how='left')['ConfirmedCases'].values
        self.train.loc[self.train.ypred0.notnull(), 'target0'] = self.train.loc[self.train.ypred0.notnull(), 'ypred0']
        self.train['ypred1'] = pd.merge(self.train[feats], yesterday[feats + ['Fatalities']], on=feats, how='left')['Fatalities'].values
        self.train.loc[self.train.ypred1.notnull(), 'target1'] = self.train.loc[self.train.ypred1.notnull(), 'ypred1']
        del self.train['ypred0'], self.train['ypred1']
        tr, vl = self.create_features(self.train.copy(), self.day)
        dvalid = xgb.DMatrix(vl[self.features0])
        ypred0 = (self.model0.predict(dvalid) + self.model1.predict(dvalid)) / 2
        dvalid = xgb.DMatrix(vl[self.features1])
        ypred1 = (self.model2.predict(dvalid) + self.model3.predict(dvalid)) / 2
        vl['ypred0'] = ypred0
        vl['ypred1'] = ypred1
        vl.loc[vl.ypred0 < vl.log0_max, 'ypred0'] = vl.loc[vl.ypred0 < vl.log0_max, 'log0_max']
        vl.loc[vl.ypred1 < vl.log1_max, 'ypred1'] = vl.loc[vl.ypred1 < vl.log1_max, 'log1_max']
        self.vl = vl
        VALID = vl[['geo', 'day', 'ypred0', 'ypred1']].copy()
        VALID.columns = ['geo', 'day', 'ConfirmedCases', 'Fatalities']
        return VALID.reset_index(drop=True)

TARGETS = ['ConfirmedCases', 'Fatalities']
def rmse(y_true, y_pred):
    return np.sqrt(mean_squared_error(y_true, y_pred))
df = pd.read_csv('../input/covid19-global-forecasting-week-4/train.csv')
df[TARGETS] = np.log1p(df[TARGETS].values)
sub_df = pd.read_csv('../input/covid19-global-forecasting-week-4/test.csv')
def preprocess(df):
    for col in ['Country_Region', 'Province_State']:
        df[col].fillna('', inplace=True)
    df['Date'] = pd.to_datetime(df['Date'])
    df['day'] = df.Date.dt.dayofyear
    df['geo'] = ['_'.join(x) for x in zip(df['Country_Region'], df['Province_State'])]
    return df
df = preprocess(df)
sub_df = preprocess(sub_df)
sub_df['day'] -= df['day'].min()
df['day'] -= df['day'].min()

TEST_FIRST = sub_df[sub_df['Date'] > df['Date'].max()]['Date'].min()
TEST_DAYS = (sub_df['Date'].max() - TEST_FIRST).days + 1
TEST_FIRST = (TEST_FIRST - df['Date'].min()).days
def get_blend(pred_dfs, weights, verbose=True):
    blend_df = pred_dfs['giba1'].copy()
    blend_df[TARGETS] = 0
    for name, pred_df in pred_dfs.items():
        blend_df[TARGETS] += weights[name] * pred_df[TARGETS].values
    return blend_df
cov_models = {'ahmet': CovidModelAhmet(), 'giba1': CovidModelGIBA(lag=1), 'giba2': CovidModelGIBA(lag=2)}
weights = {'ahmet': 0.45, 'giba1': 0.275, 'giba2': 0.275}
pred_dfs = {name: cm.predict_first_day(TEST_FIRST).sort_values('geo') for name, cm in cov_models.items()}
blend_df = get_blend(pred_dfs, weights)
eval_df = blend_df.copy()
for d in range(1, TEST_DAYS):
    pred_dfs = {name: cm.predict_next_day(blend_df).sort_values('geo') for name, cm in cov_models.items()}
    blend_df = get_blend(pred_dfs, weights)
    eval_df = eval_df.append(blend_df)

eval_df.head()