<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 We will train a model to predict drug resistance values from sequence.
This is the other general variant of supervised learning - where instead of predicting a "label" for a class (classification), we are predicting a "number" or a "value". This is called regression, analogous to, say, linear regression or logistic regression.
Step1:  Exercise
Practice what you've learned! Split the data into features, response variable, and then do a train/test split.
Step2:  Exercise
Now, let's train the Random Forest Regressor on the data.
Step3:  Make a plot of what the predictions (y-axis) against the actual values (x-axis)
Step4:  Evaluating the Model
Just as with classification tasks, we also need metrics to help evaluate how good a trained model is, given the input features.
Exercise
Look through the sklearn.metrics module. What might be a suitable metric to use?
Justify the use of two of them, and write the code that computes the evaluation metric.
Step5:  Discussion
What does the distribution of values look like? Where is its skew? How could we tell?
What would be a better way of transforming the data prior to doing ML?
Live-Coding
Let's try log10-transforming the values to be predicted.
Step6:  Note
Step7:  The r-squared score goes up. Less skew is always good.
Step8:  Challenge Exercise
Can you compare the following algorithms to see which one performs best?
RandomForestRegressor
GradientBoostingRegressor
AdaBoostRegressor
ExtraTreesRegressor
Statistical Practices
Here, I will show you how to use the ShuffleSplit iterator, alongside the cross_val_score method, to evaluate the performance of different regressor models.
The methods here aren't limited to regression, though. You can easily implement these for classification as well.
<ASSISTANT_TASK:>
Python Code:
# Load the sequence data as a Pandas dataframe.
seqids = [s.id for s in SeqIO.parse('data/hiv-protease-sequences-expanded.fasta', 'fasta')]
sequences = [s for s in SeqIO.parse('data/hiv-protease-sequences-expanded.fasta', 'fasta')]
sequences = MultipleSeqAlignment(sequences)
sequences = pd.DataFrame(np.array(sequences))
sequences.index = seqids
# Ensure that all of the letters are upper-case, otherwise the replace function in the next cell won't work.
for col in sequences.columns:
    sequences[col] = sequences[col].apply(lambda x: x.upper())
    sequences[col] = sequences[col].replace('*', np.nan)
sequences.head()
seqdf = sequences.replace(isoelectric_points.keys(), isoelectric_points.values())
seqdf.head()
# Load the drug resistance values
dr_vals = pd.read_csv('data/hiv-protease-data-expanded.csv', index_col=0)
dr_vals.set_index('seqid', inplace=True)
dr_vals.head()
# Join the sequence data together with that of one drug of interest.
drug_name = 'FPV'
data_matrix = seqdf.join(dr_vals[drug_name]).dropna()  # we have to drop NaN values because scikit-learn algorithms are not designed to accept them.
data_matrix.head()
Explanation: We will train a model to predict drug resistance values from sequence.
This is the other general variant of supervised learning - where instead of predicting a "label" for a class (classification), we are predicting a "number" or a "value". This is called regression, analogous to, say, linear regression or logistic regression.
End of explanation
# Your Answer
# Hint: to select a set of columns from a dataframe, use: dataframe[[columns]]
# Hint: the columns 0 to 98 can be expressed as a list comprehension: [i for i in range(99)]
X = data_matrix[[i for i in range(99)]]
Y = data_matrix[drug_name]
X_train, X_test, Y_train, Y_test = train_test_split(X, Y)
Explanation: Exercise
Practice what you've learned! Split the data into features, response variable, and then do a train/test split.
End of explanation
# Answer
mdl = RandomForestRegressor()
mdl.fit(X_train, Y_train)
preds = mdl.predict(X_test)
Explanation: Exercise
Now, let's train the Random Forest Regressor on the data.
End of explanation
plt.scatter(preds, Y_test)
Explanation: Make a plot of what the predictions (y-axis) against the actual values (x-axis)
End of explanation
# Metric 1: correlation coefficient.
r2_score(preds, Y_test) 
# Metric 2: mean squared error
mean_squared_error(preds, Y_test)
Explanation: Evaluating the Model
Just as with classification tasks, we also need metrics to help evaluate how good a trained model is, given the input features.
Exercise
Look through the sklearn.metrics module. What might be a suitable metric to use?
Justify the use of two of them, and write the code that computes the evaluation metric.
End of explanation
X = data_matrix[[i for i in range(99)]]
Y = data_matrix[drug_name].apply(np.log10)  # the log10 transformation is applied here.
X_train, X_test, Y_train, Y_test = train_test_split(X, Y)
mdl = RandomForestRegressor()
mdl.fit(X_train, Y_train)
preds = mdl.predict(X_test)
Explanation: Discussion
What does the distribution of values look like? Where is its skew? How could we tell?
What would be a better way of transforming the data prior to doing ML?
Live-Coding
Let's try log10-transforming the values to be predicted.
End of explanation
mean_squared_error(preds, Y_test)
Explanation: Note: the MSE goes down because of the log10 transform. However, that's exactly what we would have expected by definition.
End of explanation
r2_score(preds, Y_test) 
Explanation: The r-squared score goes up. Less skew is always good.
End of explanation
cv = ShuffleSplit(n=len(X), n_iter=10, test_size=0.3)
models = dict()
models['rf'] = RandomForestRegressor()
models['gb'] = GradientBoostingRegressor()
models['ad'] = AdaBoostRegressor()
models['ex'] = ExtraTreesRegressor()
scores = dict()
for abbr, model in models.items():
    print(abbr, model)
    score = cross_val_score(model, X, Y, cv=cv, scoring='mean_squared_error')
    scores[abbr] = -score  # a known issue in the scikit-learn package; we have to take the negative of the result.
score_summary = pd.DataFrame(scores)
score_summary = pd.DataFrame(score_summary.unstack()).reset_index()
# sns.violinplot(x=)
score_summary.columns = ['model', 'idx', 'error']
sns.violinplot(x='model', y='error', data=score_summary)
Explanation: Challenge Exercise
Can you compare the following algorithms to see which one performs best?
RandomForestRegressor
GradientBoostingRegressor
AdaBoostRegressor
ExtraTreesRegressor
Statistical Practices
Here, I will show you how to use the ShuffleSplit iterator, alongside the cross_val_score method, to evaluate the performance of different regressor models.
The methods here aren't limited to regression, though. You can easily implement these for classification as well.
End of explanation
<END_TASK>