<jupyter_start><jupyter_text>Integrate Variables & Data Preprocess<jupyter_code># =========== Pre-process ==========
df = pd.read_csv('201612-citibike-tripdata.csv')
#df = pd.read_csv('201608-citibike-tripdata.csv')
#df = all_data
# Missing values in columns
#print(df.isnull().sum())
# Snake_case the columns
def camel_to_snake(name):
    s1 = re.sub('(.)([A-Z][a-z]+)', r'\1_\2', name)
    return re.sub('([a-z0-9])([A-Z])', r'\1_\2', s1).replace(' ','').lower()
df.columns = [camel_to_snake(col) for col in df.columns]

# Parse start_time and stop_time
try:
    df['start_datetime'] = [datetime.datetime.strptime(x, '%m/%d/%Y %H:%M:%S') for x in df.start_time]
#except ValueError:
    #df['start_datetime'] = [datetime.datetime.strptime(x, '%m/%d/%Y %H:%M') for x in df.start_time]
except ValueError:
     df['start_datetime'] = [datetime.datetime.strptime(x, '%Y-%m-%d %H:%M:%S') for x in df.start_time]
df['start_day'] = [x.strftime('%Y-%m-%d') for x in df.start_datetime]

# =========== Group by date and station ===========
df['bikein'] = 1
df['bikeout'] = 1

df_bikeout = df.groupby(by=['start_stationid', 'start_day'])['bikeout'].count()
df_bikeout = pd.DataFrame(df_bikeout)
df_out = df_bikeout.reset_index(drop = False)
df_out.rename_axis({'start_stationid':'station_id'}, axis='columns', inplace=True)
print(df_out.head())

df_bikein = df.groupby(by=['end_stationid', 'start_day'])['bikein'].count()
df_bikein = pd.DataFrame(df_bikein)
df_in = df_bikein.reset_index(drop = False)
df_in.rename_axis({'end_stationid':'station_id'}, axis='columns', inplace=True)
#print(df_in.head())


# =========== Bike loss per day for each station ===========
# Merge two DFs by station_id
df_dayloss = df_out.merge(right=df_in, left_on=['station_id', 'start_day'], right_on=['station_id', 'start_day'],
                          how='outer')
df_dayloss = df_dayloss.sort_values(['station_id', 'start_day'])
df_dayloss = df_dayloss.fillna(0)
df_dayloss['bike_loss'] = df_dayloss['bikeout'] - df_dayloss['bikein']
print(df_dayloss[1:5])
df_dayloss[df_dayloss.station_id==72][1:5]
# append weekday_index column
#wekday_name=pd.to_datetime(df_dayloss.start_day).dt.weekday_name   #name
df_wekday_index=pd.to_datetime(df_dayloss.start_day).dt.weekday   # index Monday=0, Sunday=6
# Append weekday column to df_dayloss
df_dayloss['wek_index']=df_wekday_index
# =========== count docks Then append to dayloss============
jsonurl = urllib.request.urlopen('https://gbfs.citibikenyc.com/gbfs/en/station_status.json')
text = json.loads(jsonurl.read()) # <-- read from it
data=text['data']
#print(len(data['stations']))   664 stations
#print(data['stations'][0]['station_id']) # station id

# station id Extract:
total_stations_id = []
for i in range(len(data['stations'])):
    id=data['stations'][i]['station_id']
    total_stations_id.append(id)
total_stations_id=[int(x) for x in total_stations_id]
#print('station_id list:'+str(total_stations_id))

# num_bikes_available Extract:
total_num_bikes_available = []
for i in range(len(data['stations'])):
    num_bikes_available = data['stations'][i]['num_bikes_available']
    total_num_bikes_available.append(num_bikes_available)
#print('bikes availble:'+str(total_num_bikes_available))

# num_docks_available Extract:
total_docks_available = []
for i in range(len(data['stations'])):
    docks_availble=data['stations'][i]['num_docks_available']
    total_docks_available.append(docks_availble)
#print('docks availble:'+str(total_docks_available))

total_docks=[x + y for x, y in zip(total_num_bikes_available, total_docks_available)]
#print('total docks for each station:'+str(total_docks))
#to_dict

new_dict = dict(zip(total_stations_id,total_docks))
#print('zip station_id and total docks:'+str(new_dict))

# append docks to df_dayloss
df_dayloss['docks'] = df_dayloss['station_id'].map(new_dict)
# Append holiday column to df_dayloss
holiday = ['2016-01-01', '2016-01-18', '2016-02-12', '2016-02-15',
          '2016-05-08', '2016-05-30', '2016-06-19', '2016-07-04',
          '2016-09-05', '2016-10-10', '2016-11-11', '2016-11-24',
          '2016-12-26']
df_dayloss['holiday_idx'] = df_dayloss['start_day'].apply(lambda x: 1 if (x in holiday) else 0)
def weather_Data(df_climate):
    dates = []
    maxtempF = []
    cltype = []
    each_day_weather_set = []
    dk = df_climate.drop_duplicates(['start_day'])

    for x in dk['start_day']:
        dates.append(x)
        api_str = 'http://api.worldweatheronline.com/premium/v1/past-weather.ashx?key=c3e90d1cd020454594c203537172604&q=NY&date=' + x + '&format=json'       # print api_str
        time.sleep(0.02)
        r = requests.get(api_str)
        #print "request over"
        
        k = r.json()
        
        maxtempF.append(k['data']['weather'][0]['maxtempF'])

        size = len(k['data']['weather'][0]['hourly'])
        climate = []
        for l in range(size):
            climate.append(str(k['data']['weather'][0]['hourly'][l]['weatherDesc'][0].values()).split(" "))
        to_str = str(climate)
        r1=to_str.replace('[','').replace(']','').replace('"','').replace("u'",'').replace("'",'').replace(' ','')
        r2=r1.split(",")
        weather_condition_count=Counter(r2)
        most_common_weather=weather_condition_count.most_common(1)[0][0]
        each_day_weather_set.append(most_common_weather)
    new_dict = dict(zip(dates, maxtempF))
    new_dict2 = dict(zip(dates, each_day_weather_set))
    #print each_day_weather_set
    #print dates
    df_climate['High_temp'] = df_climate['start_day'].map(new_dict)
    df_climate['Climate_type'] = df_climate['start_day'].map(new_dict2)
    return df_climate
time.sleep(0.2)
weather_Data(df_dayloss)
set(df_dayloss.Climate_type)<jupyter_output><empty_output><jupyter_text>change categorical variables to dummy variables:    1. wek_index: Weekday:0; Weekend:1    2. holiday_index: Non-holiday:0; Holiday:1    3. High_temp: continuous    4. Climate_type: Snow: (1,0,0); Rain:(0,1,0); Sunny:(0,0,1); Clear:(0,0,0)<jupyter_code># Assign weekday = 0 and weekend = 1
df_dayloss['wek_index']=df_dayloss['wek_index'].replace([1,2,3,4],[0,0,0,0])
df_dayloss['wek_index']=df_dayloss['wek_index'].replace([5,6],[1,1])   # weekend
# Get Dummy Variables for 'Climate Type'
dummy_weather=pd.get_dummies(df_dayloss['Climate_type'])
df_dayloss = pd.concat([df_dayloss, dummy_weather], axis=1)

# Remove docks col and 
del df_dayloss['docks']
del df_dayloss['Climate_type']
df_dayloss[31:60]
#df_dayloss=df_dayloss.drop(df_dayloss.columns[8:11], axis=1)
# Make variables for regression in continuous columns
# cols = df_dayloss.columns.tolist()
# cols.insert(len(cols)-1, cols.pop(cols.index('Climate_type')))
# df_dayloss = df_dayloss.reindex(columns=cols)<jupyter_output><empty_output><jupyter_text>Some Visilizations    1. Average bike_loss for each station<jupyter_code>avg_loss=df_dayloss.groupby('station_id', as_index=False)['bike_loss'].mean()
plt.scatter(avg_loss.station_id, avg_loss.bike_loss,s=50)
axes = plt.gca()
axes.set_xlim([0,3500])
axes.set_ylim([-5,5])
small_sta_id=avg_loss[avg_loss.station_id<=2000]
print small_sta_id.describe()
large_sta_id=avg_loss[avg_loss.station_id>2000]
print large_sta_id.describe()<jupyter_output>station_id   bike_loss
count   291.000000  291.000000
mean    372.065292   -0.650924
std     147.551134    5.484967
min      72.000000  -49.677419
25%     290.000000   -1.903226
50%     376.000000    0.032258
75%     459.500000    1.274194
max    2000.000000   19.032258
        station_id   bike_loss
count   322.000000  322.000000
mean   3185.173913    0.605777
std     279.837081    3.272735
min    2001.000000  -11.806452
25%    3108.250000   -0.322581
50%    3238.500000    0.193548
75%    3349.750000    0.838710
max    3440.000000   33.818182<jupyter_text>Linear Regression     1. Linear regression for each station(December)    2. Coff Analysis (Larger coff represents for larger effect of that variable)    3. Significant test for each independent variable    4. R^2 & accuracy of Linear Regression Model (When test data coming in)<jupyter_code>station_id_set = set(df_dayloss['station_id'])
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    Y=a.ix[:,4]# Assign features to be X and target to be Y.
    try:
        mdl = LinearRegression().fit(X,Y)
        m = mdl.coef_
        b = mdl.intercept_
        #print 'y={}x+{}'.format(m,b)
    except:
        continue<jupyter_output><empty_output><jupyter_text>Cross Validation -- Linear Regression<jupyter_code>LR_CV_score_set = []
LR_culmilative_set = []
id_set = []
# Count RMSE (Root Mean Squared Error) for every station 
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]
    linearReg = LinearRegression()
    try:
        mean_CV_score = math.sqrt(-(cross_val_score(linearReg, X, y, cv=10, scoring='neg_mean_squared_error').mean()))
        #print 'CV score for station %s is: '%i +str(mean_CV_score)
        LR_CV_score_set.append(mean_CV_score)
        LR_culmilative_set.append(np.mean(LR_CV_score_set))
        id_set.append(i)
    except:
        continue
print 'Mean Score of Cross Validation using Linear Regression is:' + str(np.mean(LR_CV_score_set))
print 'Standard Deviation Score of Cross Validation using Linear Regression is:' + str(np.std(LR_CV_score_set))<jupyter_output>Mean Score of Cross Validation using Linear Regression is:8.2995288395
Standard Deviation Score of Cross Validation using Linear Regression is:5.35644899655<jupyter_text>Random Forest    1. Importance of Each Feature    2. Accuracy of Random Forest Model  Cross Validation --- Random Forest<jupyter_code>RF_CV_score_set = []
RF_culmilative_set = []
RF_id_set = []
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]
    rf = RandomForestRegressor(n_estimators=30,max_features=4)
    try:
        mean_CV_score = math.sqrt(-(cross_val_score(rf, X, y, cv=10, scoring='neg_mean_squared_error').mean()))
        print 'CV score for station %s is: '%i +str(mean_CV_score)
        RF_CV_score_set.append(mean_CV_score)
        RF_culmilative_set.append(np.mean(RF_CV_score_set))
        RF_id_set.append(i)
    except:
        continue
# print 'Mean Score of Cross Validation using Random Forest is: ' + str(np.mean(RF_CV_score_set))
# print 'Standard Deviation Score of Cross Validation using Random Forest is: ' + str(np.std(RF_CV_score_set))
print 'Mean Score of Cross Validation using Random Forest is: ' + str(np.mean(RF_CV_score_set))
print 'Standard Deviation Score of Cross Validation using Random Forest is: ' + str(np.std(RF_CV_score_set))<jupyter_output>Mean Score of Cross Validation using Random Forest is: 8.05785701511
Standard Deviation Score of Cross Validation using Random Forest is: 5.33366542595<jupyter_text>Predict Using Random Forest<jupyter_code>final_predict_set = []
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]
    final_model = RandomForestRegressor(n_estimators=30,max_features=4)
    final_model.fit(X,y)   # Here fit the train model 
    final_predictions = final_model.predict(X)  # I just used trained X to test. 
    final_predict_set.append(final_predictions)
final_predict_set<jupyter_output><empty_output><jupyter_text>RF better than LR !<jupyter_code># The total RMASE is counted by average rmse of all stations
# For each model(LR & RF), the culmulative RMSE changes when new station coming in.
from pylab import *
t1 = LR_culmilative_set
t2 = RF_culmilative_set
s = RF_id_set
plot(s, t1)
plot(s,t2)
xlabel('Station ID')
ylabel('Culmilative RMSE')
title('Python Line Chart: Culmulative RMSE plot for LR & RF')
grid(True)
show()<jupyter_output><empty_output><jupyter_text>Lasso   1.using alpha = 10.2<jupyter_code>from sklearn import linear_model
lasso_CV_score_set = []
lasso_culmilative_set = []
lasso_id_set = []
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]
    lasso = linear_model.Lasso(alpha=10.2)
    try:
        mean_CV_score = math.sqrt(-(cross_val_score(lasso, X, y, cv=10, scoring='neg_mean_squared_error').mean()))
        print 'CV score for station %s is: '%i +str(mean_CV_score)
        lasso_CV_score_set.append(mean_CV_score)
        lasso_culmilative_set.append(np.mean(lasso_CV_score_set))
        lasso_id_set.append(i)
    except:
        continue
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]
    lasso = linear_model.Lasso(alpha=0.1)
    lasso.fit(X,y)
for i in lasso.coef_:
    print i
#lasso.intercept_
print 'Mean Score of Cross Validation using Lasso Regression is: ' + str(np.mean(lasso_CV_score_set))
print 'Standard Deviation Score of Cross Validation using Lasso Regression is: ' + str(np.std(lasso_CV_score_set))
# grid find best hyperparameter for RF
# from sklearn.model_selection import GridSearchCV

# param_grid = [
#     {'n_estimators': [3, 10, 30], 'max_features': [2, 4, 6]},
#     {'bootstrap': [False], 'n_estimators': [3, 10], 'max_features': [2, 3, 4]},
#   ]

# forest_reg = RandomForestRegressor()

# grid_search = GridSearchCV(forest_reg, param_grid, cv=10,
#                            scoring='neg_mean_squared_error')

# box=[]
# for i in station_id_set:
#     a=df_dayloss[df_dayloss['station_id']==i]
#     X=a.ix[:,5:11]
#     y=a.ix[:,4]
#     try:
#         grid_search.fit(X, y)
#         box.append(grid_search.best_params_)
#     except:
#         continue

#Feature Importance using RF(If needed) 
from sklearn.linear_model import RandomizedLasso
from sklearn.ensemble import RandomForestRegressor
from sklearn.feature_selection import RFE
for i in station_id_set:
    a=df_dayloss[df_dayloss['station_id']==i]
    X=a.ix[:,5:16]
    y=a.ix[:,4]# Assign features to be X and target to be Y.
#     rf = RandomForestRegressor(n_estimators=30)
#     rf.fit(X, y)
#     print rf.feature_importances_
#     lr = LinearRegression()
#     rfe = RFE(lr, n_features_to_select=1)
#     rfe.fit(X,y)
    rlasso = RandomizedLasso(alpha=0.1)
    try:
        rlasso.fit(X, y)
    except:
        continue


#     rfe = RFE(rf, n_features_to_select=1)
#     rfe.fit(X,Y）
varible_names=list(df_dayloss)[5:]
# from sklearn.linear_model import RandomizedLasso
# #Feacture Importance continue
names = varible_names
names = np.array(names)
#sorted(zip(names, rf.feature_importances_),key=lambda x: x[1],reverse=True)
print 'feature imoirtance rank using lasso_score is: '
sorted(zip(names, rlasso.scores_),key=lambda x: x[1],reverse=True)
#print sorted(zip(map(lambda x: round(x, 4), rfe.ranking_), names))<jupyter_output>feature imoirtance rank using lasso_score is:<jupyter_text>training data split<jupyter_code>dayloss_72 =df_dayloss[df_dayloss['station_id'] == 79]
from sklearn.cross_validation import train_test_split
import matplotlib.pyplot as plt
%matplotlib inline
# training data split
X = dayloss_72[['wek_index','holiday_idx','docks']]
y = dayloss_72['bike_loss']
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size = 0.25, random_state = 123)
#y_train, y_test = cross_validation(y, test_size = 0.25, random_state = 123)
X_train = X_train.sort_index()
X_test = X_test.sort_index()
y_train = y_train.sort_index()
y_test = y_test.sort_index()
X_train[1:5]
plt.plot(dayloss_72.index, dayloss_72.bike_loss)<jupyter_output><empty_output><jupyter_text>Lasso<jupyter_code>from sklearn.linear_model import Lasso
from sklearn.model_selection import KFold
import numpy as np
import math
def RMSPE(y, y_pred):
    return math.sqrt( ((1- y_pred / y)**2).mean() )
def test_lasso(alphas):
    scores_lasso = []
    kf = KFold(n_splits=10, random_state = 123)
    
    for alpha in alphas:
        lasso = Lasso(alpha)
        for train, val in kf.split(X_train):
            score_cv = []
            y_pred = lasso.fit(X_train.values[train],y_train.values[train]).predict(X_train.values[val])
            score_cv.append(RMSPE(y_train.values[val], y_pred))
            score = np.mean(score_cv)
        scores_lasso.append(score)
    
    plt.figure(figsize=(14,8))
    cv_lasso = pd.Series(scores_lasso, index = alphas)
    cv_lasso.plot(title = "Validation - Lasso")
    plt.xlabel("alpha")
    plt.ylabel("RMSPE")
    print(cv_lasso.argmin())
alphas = [1e-4, 1e-3,1e-2, 0.1, 1, 5, 10, 20]
test_lasso(alphas)
alphas = [10, 12, 14, 18, 20]
test_lasso(alphas)
alphas = [10, 10.5, 11, 11.25, 11.5, 11.75, 12, 12.25, 12.5, 13]
test_lasso(alphas)
alphas = [10, 10.10,10.15, 10.20, 10.25, 10.30,10.40,10.55,10.7]
test_lasso(alphas)
lasso = Lasso(alpha=10.2)
scores_lasso = []
kf = KFold(n_splits=10, random_state = 123)
    
for alpha in alphas:
    lasso = Lasso(alpha)
    for train, val in kf.split(X_train):
        score_cv = []
        y_pred = lasso.fit(X_train.values[train],y_train.values[train]).predict(X_train.values[val])
        score_cv.append(RMSPE(y_train.values[val], y_pred))
        score = np.mean(score_cv)
    scores_lasso.append(score)
# prediction of bike loss of station 72
y_pred<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14599520944090671, 'fraction_numerical': 0.04662031898112987, 'mean_word_length': 4.723169508525577, 'pattern_counts': {'":': 0, '<': 43, '<?xml version=': 0, '>': 42, 'https://': 1, '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': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '461028', 'license': 'MIT', 'n_tokens_mistral': 7044, 'n_tokens_neox': 6517, 'n_words': 1218, 'path': '660-final_4-26.ipynb', 'repo_name': 'kkandhas/CityBike-Demand-Prediction-and-Analysis', 'size': '17117'}