<jupyter_start><jupyter_text>Read in and basic overview<jupyter_code>#path_dataset = 'gdrive/My Drive/CS109BProject/data/Capital_Projects.csv'
path_dataset = '../data/Capital_Projects.csv'

df = pd.read_csv(path_dataset)
print(df.shape)
df.replace('', np.nan, inplace=True)
df = df.dropna()
df = df[df['Budget Forecast'] > 0]
#df = df.dropna(subset=['Budget Forecast'])
df.describe(include='all')
#df.tail(200)<jupyter_output>(2259, 16)<jupyter_text>Regression Model<jupyter_code>#### (1) #### CHANGE object to datetime64[ns] for columns with dates, drop NaNs
df['Forecast Completion'] =pd.to_datetime(df['Forecast Completion'], format='%m/%d/%Y', exact=False)
df['Date Reported As Of'] = pd.to_datetime(df['Date Reported As Of'],format='%m/%d/%Y', exact=False)
df['Design Start'] = pd.to_datetime(df['Design Start'], format='%m/%d/%Y', exact=False)
date_reported_min = min(df['Date Reported As Of'])
date_min = min([min(df['Date Reported As Of']), min(df['Design Start']), min(df['Forecast Completion'])])
print(f'Minimum Date: {date_min}')
print(f'Date Reported Minimum: {date_reported_min}')
print(df.shape)
df.head()
#### (2) #### Extract MONTH and YEAR from Date Reported

df['Month'] = df['Date Reported As Of'].dt.month
df['Year'] = df['Date Reported As Of'].dt.year
month_ticks = {0:'Jan', 1:'Feb',2:'Mar',3:'Apr',4:'May',5:'Jun',6:'Jul',7:'Aug',8:'Sep',9:'Oct',10:'Nov',11:'Dec'}
fig, axes = plt.subplots(1,2, figsize=(15,6))
sns.distplot(df['Month'], ax=axes[0])
axes[0].set_xticks(list(month_ticks.keys()))
axes[0].set_xticklabels(list(month_ticks.values()))
sns.distplot(df['Year'],ax=axes[1])
plt.show()
min(df['Year'])
plt_df = df.copy()
plt_df['Forecast Completion'] = pd.to_numeric((plt_df['Forecast Completion'].copy() - min(df['Forecast Completion']))/pd.Timedelta(days=1))
plt.figure(figsize=(15,5))
sns.boxplot(
    data=plt_df,
    x='Month',
    y='Forecast Completion',
    color='red')

plt.xticks([0,1,2,3,4,5,6,7,8],['Feb', 'Mar', 'Apr','May','Jun','Jul','Sep','Oct','Dec'])
plt.show()
plt_df['Budget Forecast Norm'] = (plt_df['Budget Forecast'] - min(plt_df['Budget Forecast']))/(max(plt_df['Budget Forecast']) - min(plt_df['Budget Forecast']))
plt.figure(figsize=(15,5))
sns.boxplot(
    data=plt_df,
    x='Month',
    y='Budget Forecast Norm',
    color='red')

plt.xticks([0,1,2,3,4,5,6,7,8],['Feb', 'Mar', 'Apr','May','Jun','Jul','Sep','Oct','Dec'])
plt.show()
def get_violin_lists(df, x, colname):
  reg_data = []
  col_val = []
  min_val = min(df[x])
  for cur_val, cur_df in df.groupby(colname):
    cur_reg = pd.to_numeric((cur_df[x].copy() - min_val)/pd.Timedelta(days=1)).values.reshape(-1,1)
    reg_data.append(cur_reg)
    col_val.append(cur_val)
  return reg_data, col_val

reg_data, col_val = get_violin_lists(df, 'Forecast Completion', 'Month')
plt.figure(figsize=(18,6))
plt_vals = []
for idx in col_val:
  plt_vals.append(month_ticks[idx])
plt.xticks([0,1,2,3,4,5,6,7,8],plt_vals)
plt.ylabel('Forecast Completion')
plt.violinplot(dataset=reg_data, positions=np.arange(len(reg_data)), showmedians=True, showextrema=True)
plt.show()
reg_data, col_val = get_violin_lists(df, 'Design Start', 'Month')
plt.figure(figsize=(18,6))
plt_vals = []
for idx in col_val:
  plt_vals.append(month_ticks[idx])
plt.xticks([0,1,2,3,4,5,6,7,8],plt_vals)
plt.violinplot(dataset=reg_data, positions=np.arange(len(reg_data)), showmedians=True, showextrema=True)
plt.ylabel('Design Start')
plt.show()
def get_violin_costs(df, x, colname):
  reg_data = []
  col_val = []
  min_val = min(df[x])
  max_val = max(df[x])
  for cur_val, cur_df in df.groupby(colname):
    cur_reg = (cur_df[x]- min_val)/(max_val - min_val)
    reg_data.append(cur_reg)
    col_val.append(cur_val)
  return reg_data, col_val

reg_data, col_val = get_violin_lists(df, 'Forecast Completion', 'Month')
plt.figure(figsize=(18,6))
plt_vals = []
for idx in col_val:
  plt_vals.append(month_ticks[idx])
plt.xticks([0,1,2,3,4,5,6,7,8],plt_vals)
plt.ylabel('Budget Forecast')
plt.violinplot(dataset=reg_data, positions=np.arange(len(reg_data)), showmedians=True, showextrema=True)
plt.show()
#### (3) #### Compute stage of project using Date Reported, Final Forecast Completion
####     #### and Design Start day
### (a) Latest Schedule Changes ===>  Change Normalized to Project Duration of each project
### (b) Forecast Completion Norm ===>  Forecast Normalized to Project Duration of each project
### (c) Project Duration ===> Normalized over all projects (to normalize across big and small projects)

df['Final Date'] = df.groupby(['PID'])['Forecast Completion'].transform('last')
df['Stage'] = (0.9-0.1) * (df['Date Reported As Of'] - df['Design Start']) / (df['Final Date']- df['Design Start']) + 0.1
df['Project Duration'] = df['Final Date'] - df['Design Start']

df['Forecast Completion Norm'] = (df['Forecast Completion'] - df['Design Start'])/df['Project Duration']
df['Schedule Change Norm'] = df['Latest Schedule Changes']/(df['Project Duration']/pd.Timedelta(days=1))
max_project_duration = max(df['Project Duration'])
min_project_duration = min(df['Project Duration'])
df['Project Duration'] = (0.9 - 0.1) * (df['Project Duration'] - min_project_duration) / (max_project_duration - min_project_duration) + 0.1
print(f'Max Project Duration: {max_project_duration} \n Min Project Duration: {min_project_duration} ')
print(df.shape)
df.head()
#### (3) #### Compute stage of project using Date Reported, Final Forecast Completion
####     #### and Design Start day
### (a) Latest Budget Changes ===>  Budget Change Normalized to Final Budget of each project
### (b) Budget Forecast Norm ===> Budget Forecast Normalized to Final Budget of each project
### (c) Final Budget  ===> Normalized over all projects (to normalize across big and small projects)

df['Final Budget'] = df.groupby(['PID'])['Budget Forecast'].transform('last')
df['Budget Forecast Norm'] = df['Budget Forecast'] / df['Final Budget']
df['Budget Change Norm'] = df['Latest Budget Changes'] / df['Final Budget']
max_budget_final = max(df['Final Budget'])
min_budget_final = min(df['Final Budget'])
df['Final Budget'] = (0.9-0.1)*(df['Final Budget'] - min_budget_final)/(max_budget_final - min_budget_final) + 0.1
print(f'Max Budget: {max_budget_final} \n Min Budget: {min_budget_final}')
df.dtypes
####  (4)  #### Handle Categorical Columns using dummy variables
categorical_columns = ['Category', 'Borough', 'Current Phase'] # not including Agencies nor Decription
nyc_df = pd.get_dummies(df, columns=categorical_columns, drop_first=True)

#### (5) #### Drop redundant or useless columns, and drop NaNs
drop_columns = ['Date Reported As Of', 'PID', 'Project Name', 'Description','Managing Agency', 
                'Client Agency', 'Design Start', 'Budget Forecast', 'Latest Budget Changes', 
                'Latest Schedule Changes', 'Total Schedule Changes', 'Total Budget Changes', 
                'Forecast Completion', 'Year', 'Final Date']
nyc_df = nyc_df.drop(drop_columns, axis=1)
nyc_df.describe(include='all')

plt.figure(figsize=(15,5))
sns.boxplot(
    data=nyc_df,
    x='Month',
    y='Forecast Completion Norm',
    color='red')

plt.xticks([0,1,2,3,4,5,6,7,8],['Feb', 'Mar', 'Apr','May','Jun','Jul','Sep','Oct','Dec'])
plt.show()
max(df['Forecast Completion Norm'])
plt.figure(figsize=(15,5))
sns.boxplot(
    data=nyc_df,
    x='Month',
    y='Budget Forecast Norm',
    color='red')

plt.xticks([0,1,2,3,4,5,6,7,8],['Feb', 'Mar', 'Apr','May','Jun','Jul','Sep','Oct','Dec'])
plt.show()
nyc_df = nyc_df.dropna()
nyc_train, nyc_test = train_test_split(nyc_df, test_size=.2, random_state=90)
print(nyc_train.shape)
nyc_train.head()
qual_columns = ['Stage', 'Project Duration', 'Forecast Completion Norm', 'Schedule Change Norm',
                'Final Budget', 'Budget Forecast Norm', 'Budget Change Norm']	
#scatter_matrix(nyc_train[qual_columns].sample(frac=0.4, axis=0),alpha=0.5, figsize=(25,20));
import seaborn as sns
sns.set(style="ticks")
sns.pairplot(nyc_train[qual_columns])
nyc_train[qual_columns].corr()
f = plt.figure(figsize=(10, 7))
plt.matshow(nyc_train[qual_columns].corr(), fignum=f.number)
plt.xticks(range(len(qual_columns)), qual_columns, fontsize=14, rotation=90)
plt.yticks(range(len(qual_columns)), qual_columns, fontsize=14)
cb = plt.colorbar()
cb.ax.tick_params(labelsize=14)
plt.title('Correlation Matrix', fontsize=16);
# Budget change
y_train = nyc_train['Budget Change Norm']
X_train = nyc_train.drop(['Budget Change Norm'], axis=1)
fitted_model = LinearRegression().fit(X_train, y_train)
fitted_model.score(X_train, y_train)

# Normalization or Standardization
# DateTime => numeric - (- 2013) 
# QQ plot (normal plot)
X_train = sm.add_constant(X_train)
results = OLS(y_train, X_train).fit()
results.summary()
# Multi-collinearity
X_train_save = X_train.copy()
from statsmodels.stats.outliers_influence import variance_inflation_factor
pd.Series([variance_inflation_factor(X_train.values, i) 
               for i in range(X_train.shape[1])], 
               index=X_train.columns)
# Remove multi-collinearity
# Transforms

X_train = X_train_save.copy()
drop_mc_columns = ['Borough_Olive', 'Category_Bridges']
X_train = X_train.drop(drop_mc_columns, axis=1)
X_train['Final Budget'] = np.log(X_train['Final Budget'])
X_train['Stage'] = np.log(X_train['Stage'])
X_train['Month'] = np.log(X_train['Month'])
results = OLS(np.exp(y_train), X_train).fit()
results.summary()
# Residual Plots
y_hat_train = results.predict()
residuals = y_train - y_hat_train
fig, axes = plt.subplots(1,2,figsize=(15,5))
axes[0].set_title('Plot of Residuals')
axes[0].scatter(y_hat_train, residuals, alpha=0.1)
axes[0].axhline(0, c='black', lw=2)
axes[0].set_xlabel(r'$\hat{y}$')
axes[1].set_title('Histogram of Residual Magnitudes')
axes[1].hist(residuals, alpha=0.8)
axes[1].axhline(0, c='black', lw=2)
axes[1].set_xlabel(r'residuals')
plt.show()

import pylab 
import scipy.stats as stats

stats.probplot(residuals, dist="norm", plot=pylab)
pylab.show()<jupyter_output><empty_output><jupyter_text>GAM Models<jupyter_code>from pygam import LinearGAM, s, f
#from sklearn.metrics import r2_score
from sklearn.model_selection import KFold

smoothings = np.logspace(-3, 6, 9)
cv = 10

kf = KFold(n_splits=cv, random_state=42, shuffle=True)
scores = np.zeros((cv,len(smoothings)))

print(smoothings,'\n')
print(scores.shape, '\n')
print(nyc_train.shape)
for i, (train_index, test_index) in enumerate(kf.split(nyc_train)):
    train_df = nyc_train.iloc[train_index,:]
    test_df = nyc_train.iloc[test_index,:]
    
    for j, smoothing in enumerate(smoothings):
        cur_model = LinearGAM(s(0) + s(1) + s(2) + s(3) + s(4) + s(5) + s(6) + s(7) + s(8) + f(9) + f(10) + f(11) + f(12) + f(13) + f(14) + f(15) + f(16) + f(17) + f(18) + f(19) + f(20) + f(21) + f(22) + f(23) + f(24) + f(25) + f(26) + f(27) + f(28) + f(29) + f(30) + f(31) + f(32) + f(33) + f(34) + f(35) + f(36) + f(37) + f(38), lam = smoothing).fit(nyc_train, tra)
        scores[i,j] = r2_score(test_df.price, cur_model.predict(test_df))
                              
best_lambda = smoothings[np.argmax(np.mean(scores, axis = 0))]
best_lambda
X_train.shape, nyc_train.shape
nyc_train.dtypes
nyc_train.head()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18559556786703602, 'fraction_numerical': 0.025288178000178716, 'mean_word_length': 6.119592875318066, 'pattern_counts': {'":': 0, '<': 12, '<?xml version=': 0, '>': 20, 'https://': 0, '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': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '245270', 'license': 'MIT', 'n_tokens_mistral': 4498, 'n_tokens_neox': 4236, 'n_words': 836, 'path': 'content/project/E - Predicting Project Success/code/ProjectE_Group41_v4.ipynb', 'repo_name': 'gurdeep101/2020-CS109B', 'size': '11191'}