Untitled Python Project

#!pip install statsmodels==0.13.0

#!pip install nb_black==1.0.7

# Data Management from dateutil import relativedelta as rd import pandas as pd import numpy as np # Visualization import matplotlib.pyplot as plt import plotly as py import plotly.express as px import plotly.graph_objects as go import plotly.offline as pyo pyo.init_notebook_mode() import seaborn as sns # Regression import statsmodels.api as sm from statsmodels.formula.api import ols import statsmodels.graphics.api as smg # auto code formatting %load_ext nb_black

# Read and rename column country cty_info = pd.read_csv('country/covid19countryinfo.csv').rename(columns={'country':'Country'}) # Filter observations with aggregate country-level information # The column region for region-level observations is populated cty_info = cty_info[cty_info.region.isnull()] # Convert string data type to floating data type # Remove comma from the fields cty_info['healthexp'] = cty_info[~cty_info['healthexp'].isnull()]['healthexp'].str.replace(',','').astype('float') cty_info['gdp2019'] = cty_info[~cty_info['gdp2019'].isnull()]['gdp2019'].str.replace(',','').astype('float') # Convert to date objects with to_datetime method gov_actions = ['quarantine', 'schools', 'gathering', 'nonessential', 'publicplace'] for gov_action in gov_actions: cty_info[gov_action] = pd.to_datetime(cty_info[gov_action], format = '%m/%d/%Y') # Filter columns of interest # Note: feel free to explore other variables or datasets cty_info = cty_info[['Country', 'avghumidity', 'avgtemp', 'fertility', 'medianage', 'urbanpop', 'quarantine', 'schools', \ 'publicplace', 'gatheringlimit', 'gathering', 'nonessential', 'hospibed', 'smokers', \ 'sex0', 'sex14', 'sex25', 'sex54', 'sex64', 'sex65plus', 'sexratio', 'lung', 'femalelung', \ 'malelung', 'gdp2019', 'healthexp', 'healthperpop']] # cty_info.describe() cty_info.info() #cty_info.head(20)

<class 'pandas.core.frame.DataFrame'>
Int64Index: 175 entries, 0 to 192
Data columns (total 27 columns):
 #   Column          Non-Null Count  Dtype         
---  ------          --------------  -----         
 0   Country         175 non-null    object        
 1   avghumidity     145 non-null    float64       
 2   avgtemp         145 non-null    float64       
 3   fertility       174 non-null    float64       
 4   medianage       175 non-null    float64       
 5   urbanpop        175 non-null    float64       
 6   quarantine      24 non-null     datetime64[ns]
 7   schools         26 non-null     datetime64[ns]
 8   publicplace     14 non-null     datetime64[ns]
 9   gatheringlimit  12 non-null     float64       
 10  gathering       12 non-null     datetime64[ns]
 11  nonessential    4 non-null      datetime64[ns]
 12  hospibed        175 non-null    float64       
 13  smokers         123 non-null    float64       
 14  sex0            173 non-null    float64       
 15  sex14           173 non-null    float64       
 16  sex25           173 non-null    float64       
 17  sex54           173 non-null    float64       
 18  sex64           173 non-null    float64       
 19  sex65plus       173 non-null    float64       
 20  sexratio        172 non-null    float64       
 21  lung            163 non-null    float64       
 22  femalelung      163 non-null    float64       
 23  malelung        163 non-null    float64       
 24  gdp2019         174 non-null    float64       
 25  healthexp       170 non-null    float64       
 26  healthperpop    173 non-null    float64       
dtypes: datetime64[ns](5), float64(21), object(1)
memory usage: 38.3+ KB

# Worldometer data # ================ worldometer_data = pd.read_csv("corona-virus-report/worldometer_data.csv") # Replace missing values '' with NAN and then 0 worldometer_data = worldometer_data.replace("", np.nan).fillna(0) # Transform variables and round them up to the two decimal points # Note that there are instances of division by zero issue when there are either zero total tests or total cases worldometer_data["Case Positivity"] = round( worldometer_data["TotalCases"] / worldometer_data["TotalTests"], 2 ) worldometer_data["Case Fatality"] = round( worldometer_data["TotalDeaths"] / worldometer_data["TotalCases"], 2 ) # Resolve the division by zero issue by replacing infinity value with zero worldometer_data[worldometer_data["Case Positivity"] == np.inf] = 0 worldometer_data[worldometer_data["Case Fatality"] == np.inf] = 0 # Place case positivity into three bins worldometer_data["Case Positivity Bin"] = pd.qcut( worldometer_data["Case Positivity"], q=3, labels=["low", "medium", "high"] ) # Population Structure worldometer_pop_struc = pd.read_csv("corona-virus-report/population_structure_by_age_per_contry.csv") # Replace missing values with zeros worldometer_pop_struc = worldometer_pop_struc.fillna(0) # Merge datasets by common key country worldometer_data = worldometer_data.merge( worldometer_pop_struc, how="inner", left_on="Country/Region", right_on="Country" ) worldometer_data = worldometer_data[worldometer_data["Country/Region"] != 0] # Country information worldometer_data = worldometer_data.merge(cty_info, how="left", on="Country") # Replace space in variable names with '_' worldometer_data.columns = worldometer_data.columns.str.replace(" ", "_") worldometer_data.describe() worldometer_data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 190 entries, 0 to 189
Data columns (total 49 columns):
 #   Column                    Non-Null Count  Dtype         
---  ------                    --------------  -----         
 0   Country/Region            190 non-null    object        
 1   Continent                 190 non-null    object        
 2   Population                190 non-null    float64       
 3   TotalCases                190 non-null    int64         
 4   NewCases                  190 non-null    float64       
 5   TotalDeaths               190 non-null    float64       
 6   NewDeaths                 190 non-null    float64       
 7   TotalRecovered            190 non-null    float64       
 8   NewRecovered              190 non-null    float64       
 9   ActiveCases               190 non-null    float64       
 10  Serious,Critical          190 non-null    float64       
 11  Tot_Cases/1M_pop          190 non-null    float64       
 12  Deaths/1M_pop             190 non-null    float64       
 13  TotalTests                190 non-null    float64       
 14  Tests/1M_pop              190 non-null    float64       
 15  WHO_Region                190 non-null    object        
 16  Case_Positivity           190 non-null    float64       
 17  Case_Fatality             190 non-null    float64       
 18  Case_Positivity_Bin       190 non-null    category      
 19  Country                   190 non-null    object        
 20  Fraction_age_0-14_years   190 non-null    float64       
 21  Fraction_age_15-64_years  190 non-null    float64       
 22  Fraction_age_65+_years    190 non-null    float64       
 23  avghumidity               125 non-null    float64       
 24  avgtemp                   125 non-null    float64       
 25  fertility                 147 non-null    float64       
 26  medianage                 147 non-null    float64       
 27  urbanpop                  147 non-null    float64       
 28  quarantine                23 non-null     datetime64[ns]
 29  schools                   22 non-null     datetime64[ns]
 30  publicplace               12 non-null     datetime64[ns]
 31  gatheringlimit            12 non-null     float64       
 32  gathering                 12 non-null     datetime64[ns]
 33  nonessential              4 non-null      datetime64[ns]
 34  hospibed                  147 non-null    float64       
 35  smokers                   111 non-null    float64       
 36  sex0                      147 non-null    float64       
 37  sex14                     147 non-null    float64       
 38  sex25                     147 non-null    float64       
 39  sex54                     147 non-null    float64       
 40  sex64                     147 non-null    float64       
 41  sex65plus                 147 non-null    float64       
 42  sexratio                  146 non-null    float64       
 43  lung                      141 non-null    float64       
 44  femalelung                141 non-null    float64       
 45  malelung                  141 non-null    float64       
 46  gdp2019                   147 non-null    float64       
 47  healthexp                 146 non-null    float64       
 48  healthperpop              146 non-null    float64       
dtypes: category(1), datetime64[ns](5), float64(38), int64(1), object(4)
memory usage: 73.0+ KB

# Full data # ========= full_table = pd.read_csv("corona-virus-report/full_grouped.csv") full_table["Date"] = pd.to_datetime(full_table["Date"]) # Examine DataFrame (object type, shape, columns, dtypes) full_table.info() # type(full_table) # full_table.shape # full_table.columns # full_table.dtypes # Deep dive into the DataFrame full_table.head()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 35156 entries, 0 to 35155
Data columns (total 10 columns):
 #   Column          Non-Null Count  Dtype         
---  ------          --------------  -----         
 0   Date            35156 non-null  datetime64[ns]
 1   Country/Region  35156 non-null  object        
 2   Confirmed       35156 non-null  int64         
 3   Deaths          35156 non-null  int64         
 4   Recovered       35156 non-null  int64         
 5   Active          35156 non-null  int64         
 6   New cases       35156 non-null  int64         
 7   New deaths      35156 non-null  int64         
 8   New recovered   35156 non-null  int64         
 9   WHO Region      35156 non-null  object        
dtypes: datetime64[ns](1), int64(7), object(2)
memory usage: 2.7+ MB

# Grouped by day, country # ======================= full_grouped = pd.read_csv("corona-virus-report/full_grouped.csv") full_grouped["Date"] = pd.to_datetime(full_grouped["Date"]) # full_grouped.loc[full_grouped['Country/Region'] == 'US', 'Country/Region'] = 'USA' full_grouped.head() # Correct country names in worldometer to make them consistent with dataframe full_grouped column Country/Region before merging worldometer_data["Country/Region"].replace( { "USA": "US", "UAE": "United Arab Emirates", "S. Korea": "South Korea", "UK": "United Kingdom", }, inplace=True, ) # Draw population and country-level data full_grouped = full_grouped.merge( worldometer_data[["Country/Region", "Population"]], how="left", on="Country/Region" ) full_grouped = full_grouped.merge( cty_info, how="left", left_on="Country/Region", right_on="Country" ) full_grouped["Confirmed per 1000"] = ( full_grouped["Confirmed"] / full_grouped["Population"] * 1000 ) # Backfill data full_grouped = full_grouped.fillna(method="ffill") # Create post-invention indicators gov_actions = ["quarantine", "schools", "gathering", "nonessential", "publicplace"] for gov_action in gov_actions: full_grouped["post_" + gov_action] = ( full_grouped["Date"] >= full_grouped[gov_action] ) full_grouped["day_rel_to_" + gov_action] = ( full_grouped["Date"] - full_grouped[gov_action] ).dt.days # Create percent changes in covid19 outcomes covid_outcomes = ["Confirmed", "Deaths", "Recovered", "Active", "Confirmed per 1000"] for covid_outcome in covid_outcomes: full_grouped["pct_change_" + covid_outcome] = full_grouped.groupby( ["Country/Region"] )[covid_outcome].pct_change() full_grouped[full_grouped["pct_change_" + covid_outcome] == np.inf] = 0 # Replace space in variable names with '_' full_grouped.columns = full_grouped.columns.str.replace(" ", "_") full_grouped.info() # full_grouped.tail(20) # print(full_grouped.iloc[0,0]) # full_grouped[full_grouped['quarantine'] != None]['Country/Region'].unique() # full_grouped[full_grouped['Country/Region'] == 'Germany'][['quarantine','day_rel_to_quarantine']] # list(full_grouped.columns.values) # full_grouped.describe()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 35156 entries, 0 to 35155
Data columns (total 54 columns):
 #   Column                         Non-Null Count  Dtype  
---  ------                         --------------  -----  
 0   Date                           35156 non-null  object 
 1   Country/Region                 35156 non-null  object 
 2   Confirmed                      35156 non-null  int64  
 3   Deaths                         35156 non-null  int64  
 4   Recovered                      35156 non-null  int64  
 5   Active                         35156 non-null  int64  
 6   New_cases                      35156 non-null  int64  
 7   New_deaths                     35156 non-null  int64  
 8   New_recovered                  35156 non-null  int64  
 9   WHO_Region                     35156 non-null  object 
 10  Population                     35156 non-null  float64
 11  Country                        35156 non-null  object 
 12  avghumidity                    35156 non-null  float64
 13  avgtemp                        35156 non-null  float64
 14  fertility                      35156 non-null  float64
 15  medianage                      35156 non-null  float64
 16  urbanpop                       35156 non-null  float64
 17  quarantine                     35150 non-null  object 
 18  schools                        35124 non-null  object 
 19  publicplace                    35148 non-null  object 
 20  gatheringlimit                 35148 non-null  float64
 21  gathering                      35148 non-null  object 
 22  nonessential                   35095 non-null  object 
 23  hospibed                       35156 non-null  float64
 24  smokers                        35155 non-null  float64
 25  sex0                           35156 non-null  float64
 26  sex14                          35156 non-null  float64
 27  sex25                          35156 non-null  float64
 28  sex54                          35156 non-null  float64
 29  sex64                          35156 non-null  float64
 30  sex65plus                      35156 non-null  float64
 31  sexratio                       35156 non-null  float64
 32  lung                           35156 non-null  float64
 33  femalelung                     35156 non-null  float64
 34  malelung                       35156 non-null  float64
 35  gdp2019                        35156 non-null  float64
 36  healthexp                      35156 non-null  float64
 37  healthperpop                   35156 non-null  float64
 38  Confirmed_per_1000             35156 non-null  float64
 39  post_quarantine                35156 non-null  object 
 40  day_rel_to_quarantine          35150 non-null  float64
 41  post_schools                   35156 non-null  object 
 42  day_rel_to_schools             35124 non-null  float64
 43  post_gathering                 35156 non-null  object 
 44  day_rel_to_gathering           35148 non-null  float64
 45  post_nonessential              35156 non-null  object 
 46  day_rel_to_nonessential        35095 non-null  float64
 47  post_publicplace               35156 non-null  object 
 48  day_rel_to_publicplace         35148 non-null  float64
 49  pct_change_Confirmed           27152 non-null  float64
 50  pct_change_Deaths              21438 non-null  float64
 51  pct_change_Recovered           23768 non-null  float64
 52  pct_change_Active              25665 non-null  float64
 53  pct_change_Confirmed_per_1000  26741 non-null  float64
dtypes: float64(33), int64(7), object(14)
memory usage: 14.8+ MB

# Visualize the missingness isue in the dataset sns.heatmap(cty_info.isnull(), cbar=False)

# Create a function to plot (reusing from the previous BootCamp) def gt_n(n): countries = full_grouped[full_grouped["Confirmed"] > n]["Country/Region"].unique() temp = full_table[full_table["Country/Region"].isin(countries)] temp = temp.groupby(["Country/Region", "Date"])["Confirmed"].sum().reset_index() temp = temp[temp["Confirmed"] > n] temp["Log Confirmed"] = np.log(1 + temp["Confirmed"]) # print(temp.head()) min_date = temp.groupby("Country/Region")["Date"].min().reset_index() min_date.columns = ["Country/Region", "Min Date"] # print(min_date.head()) from_nth_case = pd.merge(temp, min_date, on="Country/Region") from_nth_case["Date"] = pd.to_datetime(from_nth_case["Date"]) from_nth_case["Min Date"] = pd.to_datetime(from_nth_case["Min Date"]) from_nth_case["N days"] = ( from_nth_case["Date"] - from_nth_case["Min Date"] ).dt.days # print(from_nth_case.head()) fig = px.line( from_nth_case, x="N days", y="Confirmed", color="Country/Region", title="N days from " + str(n) + " case", height=600, ) fig.show() fig = px.line( from_nth_case, x="N days", y="Log Confirmed", color="Country/Region", title="N days from " + str(n) + " case", height=600, ) fig.show()

def graph_cty_exceeding_cases(n, full_grouped_df, full_table_df): """ Function graphs the countries with more than n confirmed cases Parameters: n (int): The threshold minimum number of cases fully_grouped_df (pandas.DataFrame): the data source used. From Kaggle's full_grouped data set full_data_df (pandas.DataFrame): the data source used. From Kaggle's covid_19_clean_complete data set Returns: null """ # Data Preprocessing countries = full_grouped_df[full_grouped_df["Confirmed"] > n][ "Country/Region" ].unique() temp = full_table_df[full_table_df["Country/Region"].isin(countries)] temp = temp.groupby(["Country/Region", "Date"])["Confirmed"].sum().reset_index() temp = temp[temp["Confirmed"] > n] temp["Log Confirmed"] = np.log(1 + temp["Confirmed"]) # print(temp.head()) min_date = temp.groupby("Country/Region")["Date"].min().reset_index() min_date.columns = ["Country/Region", "Min Date"] # print(min_date.head()) from_nth_case = pd.merge(temp, min_date, on="Country/Region") from_nth_case["Date"] = pd.to_datetime(from_nth_case["Date"]) from_nth_case["Min Date"] = pd.to_datetime(from_nth_case["Min Date"]) from_nth_case["N days"] = ( from_nth_case["Date"] - from_nth_case["Min Date"] ).dt.days # print(from_nth_case.head()) fig = px.line( from_nth_case, x="N days", y="Confirmed", color="Country/Region", title="N days from " + str(n) + " case", height=600, ) fig.show() fig = px.line( from_nth_case, x="N days", y="Log Confirmed", color="Country/Region", title="N days from " + str(n) + " case", height=600, ) fig.show()

full_grouped.head()

Stringency_index = pd.read_csv("covid-stringency-index.csv") Stringency_index["Date"] = pd.to_datetime(Stringency_index["Day"]) df = Stringency_index.sort_values(by=['Date','Entity']) df.reset_index(inplace=True) df["Entity"].replace({"United States": "US",},inplace=True) df["Date"]=df["Date"].astype(object) df.head()

n_full_grouped = full_grouped.merge( df[["Date", "Entity", "stringency_index"]], how="inner", left_on=["Date", 'Country/Region'], right_on=["Date", "Entity"])

n_full_grouped.tail()

def plot_gov_action(covid_outcome, gov_action, full_grouped_df): """ Function plots the government action and outcome Parameters: covid_outcome (str): The outcome from covid gov_action (str): The government action to be analysed full_grouped_df (pandas.DataFrame): the data source used. From Kaggle's full_grouped data set Returns: null """ fig = px.scatter( full_grouped_df[full_grouped_df[gov_action] != None], x="day_rel_to_" + gov_action, y=covid_outcome, color="Country/Region", title="N days from " + gov_action, height=600, ) fig.update_layout(yaxis=dict(range=[0, 10])) fig.show()

gov_actions = ['quarantine', 'schools', 'gathering', 'nonessential', 'publicplace'] for i in gov_actions: plot_gov_action('pct_change_Confirmed_per_1000', i, n_full_grouped)

fig = px.scatter( n_full_grouped[n_full_grouped['stringency_index'] != None], x=n_full_grouped['stringency_index'], y='pct_change_Confirmed_per_1000', color="Country/Region", title='stringency_index', height=600, ) fig.update_layout(yaxis=dict(range=[0, 10])) fig.show()

n_full_grouped['Confirmed_per_1000'].describe()

n_full_grouped['log_Confirmed_per_1000'] = np.log(n_full_grouped['Confirmed_per_1000']+1) n_full_grouped['log_Confirmed_per_1000'].describe()

#Plot pairplot with countries organized by WHO regions g = sns.pairplot(n_full_grouped[['log_Confirmed_per_1000', 'avghumidity', 'avgtemp', 'urbanpop', 'WHO_Region']], hue='WHO_Region')

#plt.matshow(full_grouped.corr()) #plt.show() f, ax = plt.subplots(figsize=(10, 8)) corr = n_full_grouped.corr() # Generate a mask for the upper triangle mask = np.triu(np.ones_like(corr, dtype=bool)) # Generate a custom diverging colormap cmap = sns.diverging_palette(230, 20, as_cmap=True) # Plot heatmap sns.heatmap(corr, mask=mask, cmap=cmap, square=True, ax=ax)

n_full_grouped.columns

# 'post_schools', 'post_gathering', 'post_nonessential', 'post_publicplace', # Create interaction term n_full_grouped['quarXurbanpop'] = n_full_grouped['post_quarantine'] * n_full_grouped['urbanpop'] # OLS regression y = n_full_grouped['log_Confirmed_per_1000'] X = n_full_grouped[['post_quarantine', 'avghumidity', 'avgtemp', 'urbanpop', 'quarXurbanpop','stringency_index']] X = sm.add_constant(X) ols_model=sm.OLS(y,X.astype(float), missing='drop') result=ols_model.fit() print(result.summary2())

                     Results: Ordinary least squares
=========================================================================
Model:              OLS                    Adj. R-squared:     0.316     
Dependent Variable: log_Confirmed_per_1000 AIC:                46438.1512
Date:               2021-10-09 10:04       BIC:                46496.3550
No. Observations:   30177                  Log-Likelihood:     -23212.   
Df Model:           6                      F-statistic:        2329.     
Df Residuals:       30170                  Prob (F-statistic): 0.00      
R-squared:          0.317                  Scale:              0.27274   
--------------------------------------------------------------------------
                      Coef.   Std.Err.     t      P>|t|    [0.025   0.975]
--------------------------------------------------------------------------
const                 0.2098    0.0258    8.1376  0.0000   0.1593   0.2604
post_quarantine      -0.1679    0.0232   -7.2266  0.0000  -0.2134  -0.1223
avghumidity          -0.0009    0.0002   -4.3431  0.0000  -0.0014  -0.0005
avgtemp              -0.0056    0.0003  -17.4168  0.0000  -0.0062  -0.0050
urbanpop             -0.0001    0.0003   -0.4373  0.6619  -0.0007   0.0004
quarXurbanpop         0.0131    0.0003   41.1767  0.0000   0.0125   0.0137
stringency_index     -0.0017    0.0002  -11.2052  0.0000  -0.0020  -0.0014
-------------------------------------------------------------------------
Omnibus:               7653.551        Durbin-Watson:           1.448    
Prob(Omnibus):         0.000           Jarque-Bera (JB):        20774.461
Skew:                  1.356           Prob(JB):                0.000    
Kurtosis:              6.027           Condition No.:           1222     
=========================================================================
* The condition number is large (1e+03). This might indicate
strong multicollinearity or other numerical problems.