from pandas import *
import matplotlib.pyplot as plt
import datetime
from IPython.core.display import Image, HTML
from pandas.core.datetools import *
from numpy import arange,array,ones#,random,linalg
from scipy import stats
figsize(12,10)
pandas.__version__
'0.13.1'
pandas.set_option('display.height', 500)
pandas.set_option('display.max_rows', 500)
pandas.set_option('display.max_columns', 500)
pandas.set_option('display.width', 500)
pandas.set_option('display.mpl_style', False)
eq_data = read_excel('eq_disasters_by_year®ion.xlsx', 0, index_col=None, parse_dates=[1])
eq_data.drop(['id','Affected','Impacted'],axis=1,inplace=True)
#eq_data.drop(['id','Affected','Impacted', 'Homeless','Killed','Injured'],axis=1,inplace=True)
eq_data.head()
Year | Region | Occurrence | Killed | Injured | Homeless | Damage | |
---|---|---|---|---|---|---|---|
0 | 1901-01-01 | Eastern Africa | 1 | 0 | 0 | 0 | 0 |
1 | 1901-01-01 | Eastern Asia | 1 | 18 | 0 | 0 | 0 |
2 | 1902-01-01 | Central America | 2 | 2185 | 0 | 0 | 0 |
3 | 1902-01-01 | Russian Federation | 1 | 4562 | 0 | 0 | 0 |
4 | 1903-01-01 | Western Asia | 1 | 6000 | 0 | 0 | 0 |
5 rows × 7 columns
storm_data = read_excel('storm_disasters_by_year®ion.xlsx', 0, index_col=None, parse_dates=[1])
storm_data.drop(['id','Affected','Impacted'],axis=1,inplace=True)
#storm_data.drop(['id','Affected','Impacted', 'Homeless','Killed','Injured'],axis=1,inplace=True)
storm_data.head()
Year | Region | Occurrence | Killed | Injured | Homeless | Damage | |
---|---|---|---|---|---|---|---|
0 | 1900-01-01 | Northern America | 1 | 6000 | 0 | 0 | 30000 |
1 | 1902-01-01 | South-Eastern Asia | 1 | 600 | 0 | 0 | 0 |
2 | 1903-01-01 | Caribbean | 1 | 65 | 0 | 0 | 0 |
3 | 1903-01-01 | Northern America | 1 | 98 | 0 | 0 | 0 |
4 | 1904-01-01 | Southern Asia | 1 | 0 | 0 | 0 | 0 |
5 rows × 7 columns
ind = np.arange(3)
width = .35
bar(ind, eq_data[['Killed','Injured','Homeless']].sum(), width, color='r', label='Earthquakes')
bar(ind+width, storm_data[['Killed','Injured','Homeless']].sum(), width, color='b', label='Storms')
title('All Regions')
ax = gca()
ax.set_ylabel('Count', size='14')
xticks(ind+width, ['Killed','Injured','Homeless'], size='14')
legend(loc='upper left')
<matplotlib.legend.Legend at 0x6d5fd30>
ind = np.arange(3)
width = .35
region = 'Northern America'
n_america_eq = eq_data.ix[eq_data['Region'] == region]
n_america_storm = storm_data.ix[storm_data['Region'] == region]
bar(ind, n_america_eq[['Killed','Injured','Homeless']].sum(), width, color='r', label='Earthquakes')
bar(ind+width, n_america_storm[['Killed','Injured','Homeless']].sum(), width, color='b', label='Storms')
title(region)
ax = gca()
ax.set_ylabel('Count', size='14')
xticks(ind+width, ['Killed','Injured','Homeless'], size='14')
legend(loc='upper left')
<matplotlib.legend.Legend at 0x708ce10>
General earthquake statistics
eq_dmg = eq_data.ix[eq_data['Damage'] > 1000000] # Number over $1billion
print "Number of earthquake disasters over $1billion since 1900 = ", len(eq_dmg)
eq_kill = eq_data.ix[eq_data['Killed'] > 10000] # Number over 10,000
print "Number of earthquake disasters with over 10,000 dead since 1900 = ", len(eq_kill)
eq_inj = eq_data.ix[eq_data['Injured'] > 10000] # Number over 10,000
print "Number of earthquake disasters with over 10,000 injured since 1900 = ", len(eq_inj)
eq_home = eq_data.ix[eq_data['Homeless'] > 10000] # Number over 10,000
print "Number of earthquake disasters with over 10,000 homeless since 1900 = ", len(eq_home)
Number of earthquake disasters over $1billion since 1900 = 51 Number of earthquake disasters with over 10,000 dead since 1900 = 36 Number of earthquake disasters with over 10,000 injured since 1900 = 34 Number of earthquake disasters with over 10,000 homeless since 1900 = 105
General storm statistics
storm_dmg = storm_data.ix[storm_data['Damage'] > 1000000] # Number over $1billion
print "Number of storm disasters over $1billion since 1900 = ", len(storm_dmg)
storm_kill = storm_data.ix[storm_data['Killed'] > 10000] # Number over 10,000
print "Number of storm disasters with over 10,000 dead since 1900 = ", len(storm_kill)
storm_inj = storm_data.ix[storm_data['Injured'] > 10000] # Number over 10,000
print "Number of storm disasters with over 10,000 injured since 1900 = ", len(storm_inj)
storm_home = storm_data.ix[storm_data['Homeless'] > 10000] # Number over 10,000
print "Number of storm disasters with over 10,000 homeless since 1900 = ", len(storm_home)
Number of storm disasters over $1billion since 1900 = 108 Number of storm disasters with over 10,000 dead since 1900 = 15 Number of storm disasters with over 10,000 injured since 1900 = 14 Number of storm disasters with over 10,000 homeless since 1900 = 149
Regional earthquake statistics
eq_region = eq_data.groupby('Region')
eq_region_max = eq_region.max()
eq_max_occur_region = eq_region_max['Occurrence'].idxmax()
eq_max_occur_count = eq_region['Occurrence'].max().max()
print 'Region with maximum earthquake disaster occurrences:', eq_max_occur_region, '=', eq_max_occur_count
eq_max_kill_region = eq_region_max['Killed'].idxmax()
eq_max_kill_count = eq_region['Killed'].max().max()
print 'Region with maximum earthquake disaster deaths:', eq_max_kill_region, '=', eq_max_kill_count
eq_max_inj_region = eq_region_max['Injured'].idxmax()
eq_max_inj_count = eq_region['Injured'].max().max()
print 'Region with maximum earthquake disaster injuries:', eq_max_inj_region, '=', eq_max_inj_count
eq_max_home_region = eq_region_max['Homeless'].idxmax()
eq_max_home_count = eq_region['Homeless'].max().max()
print 'Region with maximum earthquake disaster homeless population:', eq_max_home_region, '=', eq_max_home_count
eq_max_dmg_region = eq_region_max['Damage'].idxmax()
eq_max_dmg_count = eq_region['Damage'].max().max()
print 'Region with maximum earthquake disaster damage ($):', eq_max_dmg_region, '=', eq_max_dmg_count
Region with maximum earthquake disaster occurrences: Eastern Asia = 14 Region with maximum earthquake disaster deaths: Eastern Asia = 242000 Region with maximum earthquake disaster injuries: Eastern Asia = 369060 Region with maximum earthquake disaster homeless population: Southern Asia = 5150000 Region with maximum earthquake disaster damage ($): Eastern Asia = 210346400
Regional storm statistics
storm_region = storm_data.groupby('Region')
storm_region_max = storm_region.max()
storm_max_occur_region = storm_region_max['Occurrence'].idxmax()
storm_max_occur_count = storm_region['Occurrence'].max().max()
print 'Region with maximum storm disaster occurrences:', storm_max_occur_region, '=', storm_max_occur_count
storm_max_kill_region = storm_region_max['Killed'].idxmax()
storm_max_kill_count = storm_region['Killed'].max().max()
print 'Region with maximum storm disaster deaths:', storm_max_kill_region, '=', storm_max_kill_count
storm_max_inj_region = storm_region_max['Injured'].idxmax()
storm_max_inj_count = storm_region['Injured'].max().max()
print 'Region with maximum storm disaster injuries:', storm_max_inj_region, '=', storm_max_inj_count
storm_max_home_region = storm_region_max['Homeless'].idxmax()
storm_max_home_count = storm_region['Homeless'].max().max()
print 'Region with maximum storm disaster homeless population:', storm_max_home_region, '=', storm_max_home_count
storm_max_dmg_region = storm_region_max['Damage'].idxmax()
storm_max_dmg_count = storm_region['Damage'].max().max()
print 'Region with maximum storm disaster damage ($):', storm_max_dmg_region, '=', storm_max_dmg_count
Region with maximum storm disaster occurrences: Western Europe = 41 Region with maximum storm disaster deaths: Southern Asia = 300317 Region with maximum storm disaster injuries: Southern Asia = 600000 Region with maximum storm disaster homeless population: Eastern Asia = 11002376 Region with maximum storm disaster damage ($): Northern America = 158230000
Group data by year
eq_yr = eq_data.groupby('Year')
eq_yr_sum = eq_yr.sum()
storm_yr = storm_data.groupby('Year')
storm_yr_sum = storm_yr.sum()
#storm_eq_yr_sum = concat([eq_yr_sum, storm_yr_sum])
storm_eq_yr_sum = concat([eq_yr_sum.rename(columns={'Occurrence': 'Occurrence (Earthquake)',
'Killed': 'Killed (Earthquake)',
'Injured': 'Injured (Earthquake)',
'Damage': 'Damage (Earthquake)',
'Homeless': 'Homeless (Earthquake)'}),
storm_yr_sum.rename(columns={'Occurrence': 'Occurrence (Storm)',
'Killed': 'Killed (Storm)',
'Injured': 'Injured (Storm)',
'Damage': 'Damage (Storm)',
'Homeless': 'Homeless (Storm)'})], axis=1)
Monotonic correlation (Spearman) between occurrence and losses
print 'Earthquake Spearman correlation = ', eq_yr_sum['Occurrence'].corr(eq_yr_sum['Damage'], method='spearman')
print 'Storm Spearman correlation = ', storm_yr_sum['Occurrence'].corr(storm_yr_sum['Damage'], method='spearman')
print 'Earthquake Spearman correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='spearman')
print 'Storm Spearman correlation from 1980 = ', storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='spearman')
Earthquake Spearman correlation = 0.871554100714 Storm Spearman correlation = 0.476538158924 Earthquake Spearman correlation from 1980 = 0.330575499317 Storm Spearman correlation from 1980 = 0.39416171794
...with rolling mean
print 'Earthquake Spearman correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(eq_yr_sum['Damage'], 3, center=True), method='spearman')
print 'Storm Spearman correlation (rolling mean) = ', rolling_mean(storm_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['Damage'], 3, center=True), method='spearman')
print 'Earthquake Spearman correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='spearman')
print 'Storm Spearman correlation from 1980 (rolling mean) = ', rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='spearman')
Earthquake Spearman correlation (rolling mean) = 0.850376849192 Storm Spearman correlation (rolling mean) = 0.386750467353 Earthquake Spearman correlation from 1980 (rolling mean) = 0.280183656033 Storm Spearman correlation from 1980 (rolling mean) = 0.462596290729
Calcuate linear correlation (Pearson) between occurrence and losses
print 'Earthquake Pearson correlation = ', eq_yr_sum['Occurrence'].corr(eq_yr_sum['Damage'], method='pearson')
print 'Storm Pearson correlation = ', storm_yr_sum['Occurrence'].corr(storm_yr_sum['Damage'], method='pearson')
print 'Earthquake Pearson correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='pearson')
print 'Storm Pearson correlation from 1980 = ', storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='pearson')
Earthquake Pearson correlation = 0.403172662868 Storm Pearson correlation = 0.217817975297 Earthquake Pearson correlation from 1980 = 0.185048723441 Storm Pearson correlation from 1980 = 0.329198804935
...with rolling mean
print 'Earthquake Pearson correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(eq_yr_sum['Damage'], 3, center=True), method='pearson')
print 'Storm Pearson correlation (rolling mean) = ', rolling_mean(storm_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['Damage'], 3, center=True), method='pearson')
print 'Earthquake Pearson correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='pearson')
print 'Storm Pearson correlation from 1980 (rolling mean) = ', rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='pearson')
Earthquake Pearson correlation (rolling mean) = 0.537345456962 Storm Pearson correlation (rolling mean) = 0.255809796698 Earthquake Pearson correlation from 1980 (rolling mean) = 0.0431994288626 Storm Pearson correlation from 1980 (rolling mean) = 0.386696332891
Calculate monotonic correlation coefficient (Spearman) between storm and earthquake indicators
print 'Occurrence Spearman correlation = ', eq_yr_sum['Occurrence'].corr(storm_yr_sum['Occurrence'], method='spearman')
print 'Damage Spearman correlation = ', eq_yr_sum['Damage'].corr(storm_yr_sum['Damage'], method='spearman')
print 'Occurrence Spearman correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], method='spearman')
print 'Damage Spearman correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Damage'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='spearman')
Occurrence Spearman correlation = 0.725262134916 Damage Spearman correlation = 0.0892185970931 Occurrence Spearman correlation from 1980 = 0.432210830826 Damage Spearman correlation from 1980 = -0.373567608862
...with rolling mean
print 'Occurrence Spearman correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['Occurrence'], 3, center=True), method='spearman')
print 'Damage Spearman correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Damage'], 3, center=True).corr(rolling_mean(storm_yr_sum['Damage'], 3, center=True), method='spearman')
print 'Occurrence Spearman correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True), method='spearman')
print 'Damage Spearman correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='spearman')
Occurrence Spearman correlation (rolling mean) = 0.868706985653 Damage Spearman correlation (rolling mean) = 0.035643804372 Occurrence Spearman correlation from 1980 (rolling mean) = 0.676578743302 Damage Spearman correlation from 1980 (rolling mean) = -0.445014662757
Calculate linear correlation coefficient (Pearson) between storm and earthquake indicators
print 'Occurrence Pearson correlation = ', eq_yr_sum['Occurrence'].corr(storm_yr_sum['Occurrence'], method='pearson')
print 'Damage Pearson correlation = ', eq_yr_sum['Damage'].corr(storm_yr_sum['Damage'], method='pearson')
print 'Occurrence Pearson correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], method='pearson')
print 'Damage Pearson correlation from 1980 = ', eq_yr_sum['1980-01-01':'2013-01-01']['Damage'].corr(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], method='pearson')
Occurrence Pearson correlation = 0.841758839515 Damage Pearson correlation = -0.636852487933 Occurrence Pearson correlation from 1980 = 0.417389682366 Damage Pearson correlation from 1980 = -0.724621687209
...with rolling mean
print 'Occurrence Pearson correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['Occurrence'], 3, center=True), method='pearson')
print 'Damage Pearson correlation (rolling mean) = ', rolling_mean(eq_yr_sum['Damage'], 3, center=True).corr(rolling_mean(storm_yr_sum['Damage'], 3, center=True), method='pearson')
print 'Occurrence Pearson correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Occurrence'], 3, center=True), method='pearson')
print 'Damage Pearson correlation from 1980 (rolling mean) = ', rolling_mean(eq_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True).corr(rolling_mean(storm_yr_sum['1980-01-01':'2013-01-01']['Damage'], 3, center=True), method='pearson')
Plot earthquake and storm indicators
fig, ax1 = subplots()
ax1.plot(eq_yr_sum.index, eq_yr_sum['Occurrence'], label='Earthquake Occurrences', color='g')
ax1.plot(storm_yr_sum.index, storm_yr_sum['Occurrence'], label='Storm Occurrences', color='m')
ax1.tick_params(labelsize=12)
ax1.set_xlabel('Year', size='16')
ax1.set_ylabel('Occurrences (EMDAT criteria)', size='16')
legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(eq_yr_sum.index, eq_yr_sum['Damage'], label='Earthquake Losses', color='g', linestyle='--')
ax2.plot(storm_yr_sum.index, storm_yr_sum['Damage'], label='Storm Losses', color='m', linestyle='--')
ax2.ticklabel_format(style='plain', axis='y')
ax2.tick_params(labelsize=12)
ax2.set_ylabel('Losses ($)', size='16')
#xlim(['1980-01-01', '2015-01-01'])
legend(loc='upper right')
<matplotlib.legend.Legend at 0x8c10610>
...with rolling mean
fig, ax1 = subplots()
ax1.plot(eq_yr_sum.index, rolling_mean(eq_yr_sum['Occurrence'], 3), label='Earthquake Occurrences', color='g')
ax1.plot(storm_yr_sum.index, rolling_mean(storm_yr_sum['Occurrence'], 3), label='Storm Occurrences', color='m')
ax1.tick_params(labelsize=12)
ax1.set_xlabel('Year', size='16')
ax1.set_ylabel('Occurrences (EMDAT criteria)', size='16')
legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(eq_yr_sum.index, rolling_mean(eq_yr_sum['Damage'], 3), label='Earthquake Losses', color='g', linestyle='--')
ax2.plot(storm_yr_sum.index, rolling_mean(storm_yr_sum['Damage'], 3), label='Storm Losses', color='m', linestyle='--')
ax2.ticklabel_format(style='plain', axis='y')
ax2.tick_params(labelsize=12)
ax2.set_ylabel('Losses ($)', size='16')
#xlim(['1980-01-01', '2015-01-01'])
legend(loc='upper right')
<matplotlib.legend.Legend at 0x8e9c810>
Difference between storm and earthquake occurrences and losses
storm_eq_occur_diff = (storm_yr_sum['Occurrence'] - eq_yr_sum['Occurrence'] )
storm_eq_dmg_diff = (storm_yr_sum['Damage'] - eq_yr_sum['Damage'] )
...with rolling mean
storm_eq_occur_diff = rolling_mean(storm_yr_sum['Occurrence'] - eq_yr_sum['Occurrence'], 3, center=True)
storm_eq_dmg_diff = rolling_mean(storm_yr_sum['Damage'] - eq_yr_sum['Damage'], 3, center=True )
Plot storm/earthquake damage differences
fig, ax1 = subplots()
ax1.plot(storm_eq_occur_diff.index, storm_eq_occur_diff, label='Occurrences', color='g')
ax1.tick_params(labelsize=12)
ax1.set_xlabel('Year', size='16')
ax1.set_ylabel('Storm - Earthquake Occurrences (EMDAT criteria)', size='16')
legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(storm_eq_occur_diff.index, storm_eq_dmg_diff, label='Damage', color='m')
ax2.ticklabel_format(style='plain', axis='y')
ax2.tick_params(labelsize=12)
ax2.set_ylabel('Storm - Earthquake Losses ($)', size='16')
xlim(['1980-01-01', '2015-01-01'])
legend(loc='upper right')
<matplotlib.legend.Legend at 0x1123d8d0>
...now with rolling mean
fig, ax1 = subplots()
ax1.plot(storm_eq_dmg_diff.index, storm_eq_occur_diff, label='Occurrences', color='g')
ax1.tick_params(labelsize=12)
ax1.set_xlabel('Year', size='16')
ax1.set_ylabel('Storm - Earthquake Occurrences (EMDAT criteria)', size='16')
legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(storm_eq_dmg_diff.index, storm_eq_dmg_diff, label='Damage ($)', color='m')
ax2.ticklabel_format(style='plain', axis='y')
ax2.tick_params(labelsize=12)
ax2.set_ylabel('Storm - Earthquake Losses ($)', size='16')
#xlim(['1980-01-01', '2015-01-01'])
legend(loc='upper right')
<matplotlib.legend.Legend at 0x112de450>
Differences between storms and earthquakes for other indicators (rolling mean only)
storm_eq_kill_diff = rolling_mean(storm_yr_sum['Killed'] - eq_yr_sum['Killed'], 3, center=True )
storm_eq_inj_diff = rolling_mean(storm_yr_sum['Injured'] - eq_yr_sum['Injured'], 3, center=True )
storm_eq_occur_diff = rolling_mean(storm_yr_sum['Occurrence'] - eq_yr_sum['Occurrence'], 3, center=True)
storm_eq_home_diff = rolling_mean(storm_yr_sum['Homeless'] - eq_yr_sum['Homeless'], 3, center=True )
years = storm_eq_occur_diff.dropna()['1980-01-01':'2015-01-01'].index.year
occur_slope, occur_int, occur_r, occur_p, occur_std = stats.linregress(years,rolling_mean(storm_eq_occur_diff.dropna()['1980-01-01':'2015-01-01'].values, 3, center=True))
occur_trend = occur_slope*years+occur_int
dmg_slope, dmg_int, dmg_r, dmg_p, dmg_std = stats.linregress(years,rolling_mean(storm_eq_dmg_diff.dropna()['1980-01-01':'2015-01-01'].values, 3, center=True))
dmg_trend = dmg_slope*years+dmg_int
kill_slope, kill_int, kill_r, kill_p, kill_std = stats.linregress(years,rolling_mean(storm_eq_kill_diff.dropna()['1980-01-01':'2015-01-01'].values, 3, center=True))
kill_trend = kill_slope*years+kill_int
inj_slope, inj_int, inj_r, inj_p, inj_std = stats.linregress(years,rolling_mean(storm_eq_inj_diff.dropna()['1980-01-01':'2015-01-01'].values, 3, center=True))
inj_trend = inj_slope*years+inj_int
home_slope, home_int, home_r, home_p, home_std = stats.linregress(years,rolling_mean(storm_eq_home_diff.dropna()['1980-01-01':'2015-01-01'].values, 3, center=True))
home_trend = home_slope*years+home_int
Plot differences between storm and earthquake for all indicators (rolling mean)
fig, ax = subplots(5, sharex=True)
ax[0].plot(storm_eq_occur_diff.index, storm_eq_occur_diff, label='Occurrence', color='k')
ax[0].plot(storm_eq_occur_diff.dropna()['1980-01-01':'2015-01-01'].index,line, linestyle='--', linewidth=0.5, color='k')
ax[0].set_title('Differences between disaster indicators (storm - earthquake)')
ax[0].set_ylabel('Occurrences')
ax[1].plot(storm_eq_dmg_diff.index, storm_eq_dmg_diff, label='Damage', color='k')
ax[1].set_ylabel('Losses ($)')
ax[2].plot(storm_eq_kill_diff.index, storm_eq_kill_diff, label='Killed', color='k')
ax[2].set_ylabel('Deaths')
ax[3].plot(storm_eq_inj_diff.index, storm_eq_inj_diff, label='Injured', color='k')
ax[3].set_ylabel('Injuries')
ax[4].plot(storm_eq_home_diff.index, storm_eq_home_diff, label='Homeless', color='k')
ax[4].set_ylabel('Homeless')
xlim(['1980-01-01', '2015-01-01'])
(722815.0, 735599.0)
r value 0.730796106505 p_value 1.3714763094e-06 standard deviation 0.183664264308
[<matplotlib.lines.Line2D at 0x12179610>]
storm_eq_occur_diff.dropna()['1980-01-01':'2015-01-01'].values
array([ -2.66666667, 9. , 13.66666667, 18.66666667, 15.66666667, 17.33333333, 13.66666667, 19.33333333, 26.33333333, 36.66666667, 31.33333333, 28.33333333, 34. , 41.66666667, 42.66666667, 37. , 37.66666667, 37.33333333, 33.66666667, 36. , 46. , 49.66666667, 38. , 33.33333333, 43. , 50.66666667, 57. , 52. , 57. , 49.66666667, 36. , 31.66666667, 19.66666667])
plot(years, storm_eq_occur_diff)
[<matplotlib.lines.Line2D at 0x12715490>]