U.S. Opiate Prescriptions/Overdoses¶

opioids = pd.read_csv('opioids.csv')
pinfo = pd.read_csv('prescriber-info.csv')
overdoses = pd.read_csv('overdoses.csv')
opioids.head()
pinfo.head()
overdoses.head()

overdoses['Deaths'].sum()

overdoses['Deaths'] = overdoses['Deaths'].str.replace(',', '')
overdoses['Deaths'] = overdoses['Deaths'].astype(int)
overdoses['Population'] = overdoses['Population'].str.replace(',', '')
overdoses['Population'] = overdoses['Population'].astype(int)
overdoses['Deaths/Population'] = (overdoses['Deaths']/overdoses['Population'])
overdoses.head()

overdoses_plot = overdoses.copy()
overdoses_plot.set_index('State', inplace=True)
overdoses_plot.head()

from matplotlib import pyplot
import matplotlib.pyplot as plt
%matplotlib inline
fig = plt.figure()
ax = overdoses_plot[['Deaths/Population']].plot(kind='bar', title ="Death per state", figsize=(15, 10), fontsize=12)
ax.set_xlabel("State", fontsize=12)
ax.set_ylabel("Deaths/Population", fontsize=12)
plt.show();

overdoses_top_5 = overdoses_plot[['Deaths/Population']].sort_values('Deaths/Population',ascending=False).iloc[0:5]
overdoses_top_5

from matplotlib import pyplot
import matplotlib.pyplot as plt
ax = overdoses_top_5[['Deaths/Population']].plot(kind='bar', title ="Death/Population per state", figsize=(10, 5), fontsize=12)
ax.set_xlabel("State", fontsize=12)
ax.set_ylabel("Deaths/Population", fontsize=12)
plt.show();

pinfo.groupby('Opioid.Prescriber').size() / pinfo.groupby('Opioid.Prescriber').size().sum()

pd.value_counts(pinfo['Opioid.Prescriber']).plot.bar()

import matplotlib.pyplot as plt
%matplotlib inline
import pandas as pd
import numpy as np
from sklearn.cross_validation import train_test_split

overdoses = pd.read_csv('overdosesnew.csv',index_col=0)
overdoses.head()

df.head()

df = pinfo[kept_variable]
df = df.drop(df.columns[[0,3]], axis=1) 
df.head()

cols = ['Specialty']
dummies = pd.concat([pd.get_dummies(df[cols][col], drop_first = True, prefix= col) for col in df[cols]], axis=1)
dummies.head()

X = dummies
X.head()

X.head()

y = df['Opioid.Prescriber']
y.head()

n_estimators = list(range(20,220,10))
max_depth = list(range(2, 22, 2)) + [None]


target = 'Opioid.Prescriber'

def rfscore2(X,y,test_size,n_estimators,max_depth):

    X_train, X_test, y_train, y_test = train_test_split(X, 
                                                        y, test_size = test_size, random_state=42) 
    rf_params = {
             'n_estimators':n_estimators,
             'max_depth':max_depth}   # parameters for grid search
    rf_gs = GridSearchCV(RandomForestClassifier(), rf_params, cv=5, verbose=1, n_jobs=-1)
    rf_gs.fit(X_train,y_train) # training the random forest with all possible parameters
    print('GridSearch results')
    print('The best parameters on the training data are:\n',rf_gs.best_params_) # printing the best parameters
    max_depth_best = rf_gs.best_params_['max_depth']      # getting the best max_depth
    n_estimators_best = rf_gs.best_params_['n_estimators']  # getting the best n_estimators
    print("best max_depth:",max_depth_best)
    print("best n_estimators:",n_estimators_best)
    best_rf_gs = RandomForestClassifier(max_depth=max_depth_best,n_estimators=n_estimators_best) # instantiate the best model
    best_rf_gs.fit(X_train,y_train)  # fitting the best model
    best_rf_score = best_rf_gs.score(X_test,y_test) 
    print ("best score is:",round(best_rf_score,2))
    print("Is the prediction smaller (S) or larger (L) than the median:\n")
    preds = best_rf_gs.predict(X_test)
    print(['S' if p == 0 else 'L' for p in best_rf_gs.predict(X_test)])
    print("")
    print("What were the probabilities of the each result above:\n")
    print("Probabilities that the number of comments is smaller than the media for each observation are:\n")
    print([('S',round(p[0],2)) if p[0] > p[1] else ('S',round(p[0],2)) for p in best_rf_gs.predict_proba(X_test)])
    print("")
    print("Confusion Matrix:\n")
    print(pd.crosstab(pd.concat([X_test,y_test],axis=1)[target], preds, rownames=['Actual Values'], colnames=['Predicted Values']))
    print('Features and their importance:\n')
    feature_importances = pd.Series(best_rf_gs.feature_importances_, index=X.columns).sort_values().tail(10)
    print(feature_importances)
    print(feature_importances.plot(kind="barh", figsize=(6,6)))
    return


rfscore2(X,y,0.3,n_estimators,max_depth)