<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 334 = 10 + 216 + 108
Step1:  ---------------------------------------------------------------------
With all the features together
---------------------------------------------------------------------
Avg RMSE for Linear Regression
Step2:  Support Vector Regression (SVR)
RBF Kernel
Step3:  ---------------------------------------------------------------------
With only Last 2 Years Crop Production
---------------------------------------------------------------------
Linear Regression
Step4:  SVR
Step5:  ---------------------------------------------------------------------
With only Satellite Data
---------------------------------------------------------------------
Avg RMSE for Linear Regression
Step6:  Support Vector Regression (SVR)
RBF Kernel
Step7:  -----------------------------------------------------
Dimension Reduction
Step8:  ------------------------------------------------------
XGBoost
Step9:  Now lets try not to Overfit.
Learning XGBoost Model on split data
Step10:  XGB with just 2 Features
Step11:  Lets check RMSE for SVR on just 2 features
Step12:  SVR (Polynomial)
All features
Step13:  Without Satellite Features
Step14:  Ridge, Lasso on All features
For X_all
Step15:  Avg RMSE is is 700
For X_crp
Step16:  XGBoost After Dimension Reduction
<ASSISTANT_TASK:>
Python Code:
select = colss[8:226]
X = rice[select]
y = rice["Value"]*1000
X.describe()
# Z-Score Normalization 
colms = list(X.columns)
for col in colms:
    col_zscore = col + '_zscore'
    X[col_zscore] = (X[col] - X[col].mean())/X[col].std(ddof=0)
cols = list(X.columns.values)
len(cols)
# Contains all the features (Last 2 years Crop Production and the Satellite Data)
select = cols[218:436]
X_all = X[select]
# Just the last 2 years Crop Production
select1 = cols[218:220]
X_crp = X[select1]
# Just the Satellite Data
select2 = cols[220:436]
X_sat = X[select2]
Explanation: 334 = 10 + 216 + 108
End of explanation
clf = LinearRegression()
scores = cross_val_score(clf, X_all, y, cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
    
print(scores)
avg_rmse = scores.mean()
print("\nAvg RMSE is "+str(scores.mean()))
Explanation: ---------------------------------------------------------------------
With all the features together
---------------------------------------------------------------------
Avg RMSE for Linear Regression
End of explanation
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=8000.0, epsilon=0.1, kernel='rbf', gamma=0.001)
scores = cross_val_score(clf, X_all, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=2000.0, epsilon=0.1, kernel='rbf', gamma=0.008)
scores = cross_val_score(clf, X_all, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=1000000.0, epsilon=0.1, kernel='rbf', gamma=0.000001)
scores = cross_val_score(clf, X_all, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: Support Vector Regression (SVR)
RBF Kernel
End of explanation
clf = LinearRegression()
scores = cross_val_score(clf, X_crp, y, cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
    
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: ---------------------------------------------------------------------
With only Last 2 Years Crop Production
---------------------------------------------------------------------
Linear Regression
End of explanation
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=500000.0, epsilon=0.1, kernel='rbf', gamma=0.0008)
scores = cross_val_score(clf, X_crp, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=8000.0, epsilon=0.1, kernel='rbf', gamma=0.01)
scores = cross_val_score(clf, X_crp, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: SVR
End of explanation
clf = LinearRegression()
scores = cross_val_score(clf, X_sat, y, cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
    
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: ---------------------------------------------------------------------
With only Satellite Data
---------------------------------------------------------------------
Avg RMSE for Linear Regression
End of explanation
# 5 Fold CV, to calculate avg RMSE
clf = SVR(C=1000.0, epsilon=0.1, kernel='rbf', gamma=0.01)
scores = cross_val_score(clf, X_sat, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: Support Vector Regression (SVR)
RBF Kernel
End of explanation
from sklearn.decomposition import PCA
def fPCA(x): 
    pca = PCA(n_components=x)
    pcax = pca.fit_transform(X_all)
    # 5 Fold CV, to calculate avg RMSE
    clf = SVR(C=2000.0, epsilon=0.1, kernel='rbf', gamma=0.008)
    scores = cross_val_score(clf, pcax, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
    for i in range(0,5):
        scores[i] = sqrt(-1*scores[i])
    print(scores)
    avg_rmse = scores.mean()
    print("\n\nAvg RMSE is ",scores.mean())
fPCA(42)
Explanation: -----------------------------------------------------
Dimension Reduction
End of explanation
import xgboost as xgb
# prepare dict of params for xgboost to run with
xgb_params = {
    'n_trees': 500, 
    'eta': 0.005,
    'max_depth': 4,
    'subsample': 0.95,
    'objective': 'reg:linear',
    'eval_metric': 'rmse',
    'base_score': np.mean(y), # base prediction = mean(target)
    'silent': 1
}
# form DMatrices for Xgboost training
dtrain = xgb.DMatrix(X_all, y)
# xgboost, cross-validation
cv_result = xgb.cv(xgb_params, 
                   dtrain, 
                   num_boost_round=900, # increase to have better results (~700)
                   early_stopping_rounds=50,
                   verbose_eval=50, 
                   show_stdv=False
                  )
num_boost_rounds = len(cv_result)
print(num_boost_rounds)
# train model
model = xgb.train(dict(xgb_params, silent=0), dtrain, num_boost_round=num_boost_rounds)
from sklearn.metrics import r2_score
print(r2_score(dtrain.get_label(), model.predict(dtrain)))
X_train, X_test, y_train, y_test = train_test_split(X_all, y, test_size=0.2, random_state=1)
dtest = xgb.DMatrix(X_test)
y_predict = model.predict(dtest)
rmse = sqrt(mean_squared_error(y_predict, y_test))
print(rmse)
Explanation: ------------------------------------------------------
XGBoost
End of explanation
X_train, X_test, y_train, y_test = train_test_split(X_all, y, test_size=0.2, random_state=1936)
dtrain = xgb.DMatrix(X_train, y_train)
# xgboost, cross-validation
cv_result = xgb.cv(xgb_params, 
                   dtrain, 
                   num_boost_round=700, # increase to have better results (~700)
                   early_stopping_rounds=50,
                   verbose_eval=50, 
                   show_stdv=False
                  )
num_boost_rounds = len(cv_result)
print(num_boost_rounds)
# train model
model = xgb.train(dict(xgb_params, silent=0), dtrain, num_boost_round=num_boost_rounds)
print(r2_score(dtrain.get_label(), model.predict(dtrain)))
dtest = xgb.DMatrix(X_test)
y_predict = model.predict(dtest)
rmse = sqrt(mean_squared_error(y_predict, y_test))
print(rmse)
avg_rmse = (481.67 + 542.952 + 577.985 + 544.2265 + 589.76 + 494.5)/6
avg_rmse
yt = y_test.as_matrix()
p = pd.DataFrame()
p["y_predicted"] = y_predict
p["y_test"] = yt
p["y_predicted"] = p["y_predicted"].round(decimals=1)
p["y_test"] = p["y_test"].round(decimals=1)
print p
Explanation: Now lets try not to Overfit.
Learning XGBoost Model on split data
End of explanation
X_train, X_test, y_train, y_test = train_test_split(X_crp, y, test_size=0.2, random_state=1)
dtrain = xgb.DMatrix(X_train, y_train)
# xgboost, cross-validation
cv_result = xgb.cv(xgb_params, 
                   dtrain, 
                   num_boost_round=700, # increase to have better results (~700)
                   early_stopping_rounds=50,
                   verbose_eval=50, 
                   show_stdv=False
                  )
num_boost_rounds = len(cv_result)
print(num_boost_rounds)
# train model
model = xgb.train(dict(xgb_params, silent=0), dtrain, num_boost_round=num_boost_rounds)
print(r2_score(dtrain.get_label(), model.predict(dtrain)))
dtest = xgb.DMatrix(X_test)
y_predict2 = model.predict(dtest)
rmse = sqrt(mean_squared_error(y_predict2, y_test))
print(rmse)
Explanation: XGB with just 2 Features
End of explanation
X_train, X_test, y_train, y_test = train_test_split(X_crp, y, test_size=0.2, random_state=1)
clf = SVR(C=8000.0, epsilon=0.1, kernel='rbf', gamma=0.01)
clf.fit(X_train, y_train)
y_predict = clf.predict(X_test)
rmse = sqrt(mean_squared_error(y_predict, y_test))
print(rmse)
yt = y_test.as_matrix()
p = pd.DataFrame()
p["y_predicted"] = y_predict
p["y_test"] = yt
p["y_predicted"] = p["y_predicted"].round(decimals=1)
p["y_test"] = p["y_test"].round(decimals=1)
p.describe()
# print (p)
Explanation: Lets check RMSE for SVR on just 2 features
End of explanation
# 5 Fold CV, to calculate avg RMSE
clf = SVR(kernel='poly', gamma='auto', degree=4, coef0=6)
scores = cross_val_score(clf, X_all, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
    
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: SVR (Polynomial)
All features
End of explanation
# 5 Fold CV, to calculate avg RMSE
clf = SVR(kernel='poly', gamma='auto', degree=3, coef0=6)
scores = cross_val_score(clf, X_crp, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
for i in range(0,5):
    scores[i] = sqrt(-1*scores[i])
    
print(scores)
avg_rmse = scores.mean()
print("\n\nAvg RMSE is ",scores.mean())
Explanation: Without Satellite Features
End of explanation
from sklearn import linear_model
reg = linear_model.RidgeCV(alphas=[1,2,3,4,5,6,7,7.1,7.2,7.3,8,9,100,120,125,130,140,150])
reg.fit(X_all, y)
reg.alpha_
X_train, X_test, y_train, y_test = train_test_split(X_all, y, test_size=0.2, random_state=50)
reg = linear_model.Ridge(alpha = 130)
reg.fit (X_train, y_train)
y_pred = reg.predict(X_test)
rmse = sqrt(mean_squared_error(y_pred, y_test))
print(rmse)
Explanation: Ridge, Lasso on All features
For X_all
End of explanation
from sklearn import linear_model
reg = linear_model.RidgeCV(alphas=[0.1,0.9,1,1.5,2,3,4,5,6,7,7.1,7.2,7.3,8,9,100,120,125,130,140,150])
reg.fit(X_crp, y)
reg.alpha_
X_train, X_test, y_train, y_test = train_test_split(X_crp, y, test_size=0.2, random_state=50)
reg = linear_model.Ridge(alpha = 1)
reg.fit (X_train, y_train)
y_pred = reg.predict(X_test)
rmse = sqrt(mean_squared_error(y_pred, y_test))
print(rmse)
# data to plot
n_groups = 4
rmse_crp = (800, 750, 750, 750)
rmse_sat = (700, 830, 830, 540)
 
# create plot
# fig, ax = plt.subplots()
index = np.arange(n_groups)
bar_width = 0.35
opacity = 0.8
 
plt.figure(figsize=(15,10))
rects1 = plt.bar(index, rmse_crp, bar_width,
                 alpha=opacity,
                 color='b',
                 label='Without Satellite Data')
 
rects2 = plt.bar(index + bar_width, rmse_sat, bar_width,
                 alpha=opacity,
                 color='g',
                 label='With Satellite Data')
 
plt.xlabel('ML Algorithms', fontsize=35)
plt.ylabel('RMSE', fontsize=35)
# plt.title('Effect of Satellite Data', fontsize=45)
plt.xticks(index + bar_width/2, ('Linear Regression \n (Ridge, Lasso)', 'SVR (Gaussian Kernel)', 'SVR (Polynomial)', 'XGBoost'), fontsize=17)
plt.legend(fontsize=15)
 
# plt.tight_layout()
plt.show()
Explanation: Avg RMSE is is 700
For X_crp
End of explanation
pca = PCA(n_components=42)
pcax = pca.fit_transform(X_all)
X_train, X_test, y_train, y_test = train_test_split(pcax, y, test_size=0.2)
dtrain = xgb.DMatrix(X_train, y_train)
# xgboost, cross-validation
cv_result = xgb.cv(xgb_params, 
                   dtrain, 
                   num_boost_round=700, # increase to have better results (~700)
                   early_stopping_rounds=50,
                   verbose_eval=50, 
                   show_stdv=False
                  )
num_boost_rounds = len(cv_result)
print(num_boost_rounds)
# train model
model = xgb.train(dict(xgb_params, silent=0), dtrain, num_boost_round=num_boost_rounds)
print(r2_score(dtrain.get_label(), model.predict(dtrain)))
dtest = xgb.DMatrix(X_test)
y_predict = model.predict(dtest)
rmse = sqrt(mean_squared_error(y_predict, y_test))
print(rmse)
# # y_pred = model.predict(dtest)
# # 5 Fold CV, to calculate avg RMSE
# # clf = SVR(C=1000.0, epsilon=0.1, kernel='rbf', gamma=0.01)
# scores = cross_val_score(model, X_all, y.values.ravel(), cv=5, scoring='neg_mean_squared_error')
# for i in range(0,5):
#     scores[i] = sqrt(-1*scores[i])
    
# print(scores)
# avg_rmse = scores.mean()
# print("\n\nAvg RMSE is ",scores.mean())
Explanation: XGBoost After Dimension Reduction
End of explanation
<END_TASK>