Answer all questions and submit them either as an IPython notebook, LaTeX document, or Markdown document. Provide full answers for each question, including interpretation of the results. Each question is worth 25 points.

The titanic.xls spreadsheet in the data directory contains data regarding the passengers on the Titanic when it sank in 1912. A recent Kaggle competition was based on predicting survival for passengers based on the attributes in the passenger list.

Use scikit-learn to build both a support vector classifier and a logistic regression model to predict survival on the Titanic. Use cross-validation to assess your models, and try to tune them to improve performance.

Discuss the benefits and drawbacks of both approaches for application to such problems.

parameters2={'penalty':('l1','l2'),'C':[0.0001,0.001,0.01,0.1,1,10,100,1000,10000],'intercept_scaling':[1,10,50,100,500,1000]}# increase intercept_scaling to minimize the penalty on the interceptnp.random.seed(3535)gscv2=GridSearchCV(LogisticRegression(),parameters2,cv=5)# default is 3-foldgscv2.fit(X_normalized,y)

The file TNNASHVI.txt in your data directory contains daily temperature readings for Nashville, courtesy of the Average Daily Temperature Archive. This data, as one would expect, oscillates annually. Using either GPy or PyMC, use a Gaussian process to fit a non-parametric regression model to this data, choosing an appropriate covariance function. Plot 10 regression lines drawn from your process.

/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/matplotlib/figure.py:1742: UserWarning:This figure includes Axes that are not compatible with tight_layout, so its results might be incorrect.

In [184]:

draws=m.posterior_samples(X,full_cov=True,size=10)

/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/GPy/core/gp.py:492: RuntimeWarning:covariance is not positive-semidefinite.

The data in prostate.data.txt come from a study by Stamey et al. (1989), which examined the correlation between the level of prostate-specific antigen (lpsa) and a number of clinical measures in men who were about to receive a radical prostatectomy. The variables are log cancer volume (lcavol), log prostate weight (lweight), age, log of the amount of benign prostatic hyperplasia (lbph), seminal vesicle invasion (svi), log of capsular penetration (lcp), Gleason score (gleason), and percent of Gleason scores 4 or 5 (pgg45).

An alternative method for model assessment is to fit the models on a set of bootstrap samples, and then keep track of how well it predicts the original training set. If $\hat{f}^b(x_i)$ is the predicted value at $x_i$, from the model fitted to the bth bootstrap dataset, such an estimate is:
$$\frac{1}{B} \frac{1}{N} \sum_{b=1}^B \sum_{i=1}^N L(y_i,\hat{f}^b(x_i)) $$
However, because the bootstrap samples tend to contain many observations in common among the set of bootstrap samples, this estimate will tend to underestimate the true error rate. The so-called .632 estimator aleviates this bias by returning a weighted average of the training error (average loss over the training sample) and the leave-one-out (LOO) bootstrap error:
$$\hat{err}^{(.632)} = 0.368 \, \bar{err} + 0.632 \, \hat{err}^{(1)}$$
where:
$$\bar{err} = \frac{1}{N}\sum_{i=1}^N L(y_i, \hat{f}(x_i)) $$
Repeat the assesment from part (1) using the .632 estimator, and compare the result to the other approaches.

In [46]:

# Part 1prostate=pd.read_table('/Users/fonnescj/Teaching/Bios8366/data/prostate.data.txt',index_col=0)prostate.head()

# Part 2B=50scores_632=[]fori,xvarsinenumerate(subsets):# Calculate bootstrap sample indicesX=prostate[xvars].valuesy_=y.valuesbootstrap_ind=[np.random.choice(range(len(X)),len(X))for_inrange(B)]# Fit models to bootstrap samplesregmod=linear_model.LinearRegression()bootstrap_fits=[regmod.fit(X[i],y_[i])foriinbootstrap_ind]# Calculate loss for each observation over models that did not use it to fit the modelloss=[np.mean([(regmod.predict(X[i].reshape(1,-1))-y_[i])**2forbinrange(B)ifinotinbootstrap_ind[b]])foriinrange(len(X))]# LOO errorerr_1=np.mean(loss)# Naive (training) errorerr_bar=((regmod.fit(X,y_).predict(X)-y_)**2).mean()scores_632.append(err_bar)

Fit a series of random-forest classifiers to the very low birthweight infant data (vlbw.csv), to explore the sensitivity to the parameter m, the number of variables considered for splitting at each step. Plot both the out-of-bag error as well as the test error against a suitably-chosen range of values for m.

/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "
/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "
/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "
/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "
/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "
/Users/fonnescj/anaconda3/envs/bios8366/lib/python3.5/site-packages/sklearn/ensemble/forest.py:439: UserWarning: Some inputs do not have OOB scores. This probably means too few trees were used to compute any reliable oob estimates.
warn("Some inputs do not have OOB scores. "

Use a grid search to optimize the number of estimators and max_depth for a Gradient Boosted Decision tree using the very low birthweight infant data. Plug this optimal max_depth into a single decision tree. Does this single tree over-fit or under-fit the data? Repeat this for the Random Forest. Construct a single decision tree using the max_depth which is optimal for the Random Forest. Does this single tree over-fit or under-fit the data?