<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Sale price distribution
First step is to look at the target sale price for the training data set, i.e. the column we're trying to predict.
Step1:  The sale price is in hte hundreds of thousands, so let's divide the price by 1000 to get more manageable numbers.
Step2:  The distribution is skewed (as demonstrated by the large z-score (and small pvalue) of teh skewtest). It is right skewed (the skew is positive). Skewed distribution are not ideal for linear models, which often assume a normal distribution. One way to correct for right-skewness is to take the log [1,2,3]
[1] http
Step3:  Merge the training and test datasets for data preparation
We're going to explore the training dataset and apply some transformations to it (fixing missing values, transforming columns etc). We'll need to apply the same transformations to the test dataset. To make that easy, let's put the training and test datasets into one dataframe.
Step4:  Features
The dataset is wide with 78 features.
Step5:  We've got 3 data types
Step6:  Split the data between categorical and numerical features
Step7:  Numerical features
Create a numerical dataset to keep track of the features
Step8:  We've got 36 numerical features. We can use the describe method to get some statistics
Step9:  But that's a lot of numbers to digest. Better get started plotting! To help with plotting, but also to improve linear regression models, we're going to standardize our data. But before that we must deal with the NaN values.
http
Step10:  Based on the description, the null values for the MasVnrArea should be 0 (no massonry veneer type)
Step11:  For the GarageYrBlt, replace by the year the house was built.
Step12:  Standardize the data
Step14:  Plot violinplots for each feature
The violin plots give us some idea of the distribution of data for each feature. We can look for things like skewness, non-normality, and the presence of outliers.
Step15:  Many of the features are higly skewed with very long tails.
Step16:  Most of these are right skewed as well. BsmtFullBath has some discrete values (number of bathrooms).
Step17:  Some features, such as BsmtFinSF2, are almost constant (blobs with long tail) as can be seen below
Step18:  Drop nearly constant features
Step19:  Log transform the other features if they have a high skewness
Using a log transformation for some of the skewed features should help, as illustrated below. We use the raw data (not the standardized one) because we need positive values for the log function (we'll standardize the transformed variables later).
Step20:  Check the sign of the skewness for all these
Step22:  Let's apply a log1p transform to all these and plot the distributions again
Step23:  Now our originally skewed features look more symmetric.
Step24:  Save transformed numerical data
Use the storage magic to communicate between notebooks.
Step25:  Feature selection
We're now in a good position to identify the key numerical features. Those should be hightly correlated with the sale price.
Step26:  Let's keep only the features that have a high enough correlation with the price (correlation less than 0.2)
<ASSISTANT_TASK:>
Python Code:
target = pd.read_csv('../data/train_target.csv')
target.describe()
Explanation: Sale price distribution
First step is to look at the target sale price for the training data set, i.e. the column we're trying to predict.
End of explanation
target = target / 1000
sns.distplot(target);
plt.title('SalePrice')
import scipy as sp
sp.stats.skew(target)
sp.stats.skewtest(target)
Explanation: The sale price is in hte hundreds of thousands, so let's divide the price by 1000 to get more manageable numbers.
End of explanation
logtarget = np.log1p(target)
print('skewness of logtarget = ', sp.stats.skew(logtarget)[0])
print('skewness test of logtarget = ', sp.stats.skewtest(logtarget))
sns.distplot(logtarget)
plt.title(r'log(1 + SalePrice)')
Explanation: The distribution is skewed (as demonstrated by the large z-score (and small pvalue) of teh skewtest). It is right skewed (the skew is positive). Skewed distribution are not ideal for linear models, which often assume a normal distribution. One way to correct for right-skewness is to take the log [1,2,3]
[1] http://fmwww.bc.edu/repec/bocode/t/transint.html 
[2] https://www.r-statistics.com/2013/05/log-transformations-for-skewed-and-wide-distributions-from-practical-data-science-with-r/
[3] Alexandru Papiu's notebook https://www.kaggle.com/apapiu/house-prices-advanced-regression-techniques/regularized-linear-models/commentsnotebook 
We apply the function $x \rightarrow \log(1 + x)$ because it is always positive for $x \geq 0$
End of explanation
raw_train = pd.read_csv('../data/train_prepared_light.csv')
raw_test = pd.read_csv('../data/test_prepared_light.csv')
df = pd.concat([raw_train, raw_test], keys=['train', 'test'])
df.shape
ncategories = sum(df.dtypes == object)
ncategories
df.head()
df.tail()
Explanation: Merge the training and test datasets for data preparation
We're going to explore the training dataset and apply some transformations to it (fixing missing values, transforming columns etc). We'll need to apply the same transformations to the test dataset. To make that easy, let's put the training and test datasets into one dataframe.
End of explanation
df.columns, len(df.columns)
Explanation: Features
The dataset is wide with 78 features.
End of explanation
df.dtypes.unique()
Explanation: We've got 3 data types: int, float and object
End of explanation
is_categorical = (df.dtypes == object)
is_numerical = ~is_categorical
Explanation: Split the data between categorical and numerical features
End of explanation
dfnum = df.loc[:, is_numerical].copy()
dfnum.columns, len(dfnum.columns)
Explanation: Numerical features
Create a numerical dataset to keep track of the features
End of explanation
dfnum.describe()
Explanation: We've got 36 numerical features. We can use the describe method to get some statistics:
End of explanation
cols_with_nulls = dfnum.columns[dfnum.isnull().sum() > 0]
cols_with_nulls
dfnum.shape
dfnum[cols_with_nulls].isnull().sum().sort_values(ascending=False)
#.plot(kind='bar')
Explanation: But that's a lot of numbers to digest. Better get started plotting! To help with plotting, but also to improve linear regression models, we're going to standardize our data. But before that we must deal with the NaN values.
http://sebastianraschka.com/Articles/2014_about_feature_scaling.html
Deal with NaN values
End of explanation
# We may want to refine this in the future. Perhaps build a model to predict the missing GarageCars from the other features?
median_list = 'LotFrontage', 'BsmtFullBath','BsmtHalfBath', 'GarageCars', 'GarageArea'
zero_list = 'MasVnrArea', 'BsmtFinSF1', 'BsmtFinSF2', 'TotalBsmtSF', 'BsmtUnfSF'
for feature in median_list:
    dfnum[feature].fillna(dfnum[feature].median(), inplace=True) 
for feature in zero_list:
    dfnum[feature].fillna(0, inplace=True)
Explanation: Based on the description, the null values for the MasVnrArea should be 0 (no massonry veneer type)
End of explanation
dfnum.GarageYrBlt.fillna(dfnum.YearBuilt[dfnum.GarageYrBlt.isnull()], inplace=True)
# Check that we got rid of the nulls
dfnum.isnull().sum().any()
# Assign to the slice (see the copy / write problem in Pandas)
df.loc[:, is_numerical] = dfnum
Explanation: For the GarageYrBlt, replace by the year the house was built.
End of explanation
def standardize(df):
    return sk.preprocessing.StandardScaler().fit_transform(df)
dfnum_t = dfnum.apply(standardize)
dfnum_t.head()
Explanation: Standardize the data
End of explanation
def violinplot(df, ax=None):
    if ax is None:
        ax = plt.gca()
    sns.violinplot(df, ax=ax)
    for xlab in ax.get_xticklabels():
        xlab.set_rotation(30)
        
def featureplot(df, nrows=1, figsize=(12,8), plotfunc=violinplot):
    Plot the dataframe features
    width, height = figsize
    fig, axes = plt.subplots(nrows, 1, figsize=(width, height * nrows));
    i = 0
    plots_per_figure = df.shape[1] // nrows
    if nrows == 1:
        axes = [axes]
    for j, ax in zip(range(plots_per_figure, df.shape[1] + 1, plots_per_figure), axes):
        plotfunc(df.iloc[:, i:j], ax=ax)
        i = j
dfnum_t.head()
train = dfnum.loc['train',:]
train_t = dfnum_t.loc['train',:]
Explanation: Plot violinplots for each feature
The violin plots give us some idea of the distribution of data for each feature. We can look for things like skewness, non-normality, and the presence of outliers.
End of explanation
featureplot(train_t.iloc[:, 0:9])
Explanation: Many of the features are higly skewed with very long tails.
End of explanation
featureplot(train_t.iloc[:, 9:18])
Explanation: Most of these are right skewed as well. BsmtFullBath has some discrete values (number of bathrooms).
End of explanation
fig, ax = plt.subplots(1,1, figsize=(4, 4))
sns.distplot(train_t['BsmtFinSF2'], ax=ax)
ax.set_title('Distribution of BsmtFinSF2')
Explanation: Some features, such as BsmtFinSF2, are almost constant (blobs with long tail) as can be seen below
End of explanation
def test_nearly_constant(series):
    counts = series.value_counts()
    max_val_count = max(counts)
    other_val_count = counts.drop(counts.argmax()).sum()
    return other_val_count / max_val_count < 0.25
is_nearly_constant = train_t.apply(test_nearly_constant)
is_nearly_constant.value_counts()
dropme = train_t.columns[is_nearly_constant]
dropme
df = df.drop(dropme, axis=1)
train = train.drop(dropme, axis=1)
train_t = train_t.drop(dropme, axis=1)
Explanation: Drop nearly constant features
End of explanation
fig, axes = plt.subplots(1,2, figsize=(8, 4))
sns.distplot(train['LotArea'], ax=axes[0])
sns.distplot(np.log1p(train['LotArea']), ax=axes[1])
zfactors = sp.stats.skewtest(train)[0]
sns.distplot(zfactors)
is_skewed = np.abs(zfactors) > 10
pd.Series(data=zfactors, index=train.columns)[is_skewed].sort_values().plot(kind='barh')
plt.title('Z-factor for skewtest')
Explanation: Log transform the other features if they have a high skewness
Using a log transformation for some of the skewed features should help, as illustrated below. We use the raw data (not the standardized one) because we need positive values for the log function (we'll standardize the transformed variables later).
End of explanation
assert all(np.sign(sp.stats.skew(train)[is_skewed]) > 0)
Explanation: Check the sign of the skewness for all these
End of explanation
def transform_skewed_colums(dfnum, is_skewed=is_skewed):
    
    dfnum: dataframe to transform
    dropme: columns to drop
    is_skewed: iterable of length dfnum.columns indicating if a column is skewed
    
    dfnum2 = dfnum.copy()
    for feature, skewed_feature in zip(dfnum.columns, is_skewed):
        if skewed_feature:
            dfnum2[feature] = np.log1p(dfnum[feature])
    dfnum_t2 = dfnum2.apply(standardize)
    return dfnum_t2
# the transformed dataset has fewer columns and we only want those
dfnum_t2 = transform_skewed_colums(df.loc[:, is_numerical])
dfnum_t2.iloc[:, is_skewed].columns
zfactors2 = sp.stats.skewtest(dfnum_t2)[0]
pd.Series(data=zfactors2, index=dfnum_t2.columns)[is_skewed].sort_values().plot(kind='barh')
Explanation: Let's apply a log1p transform to all these and plot the distributions again
End of explanation
featureplot(dfnum_t2.iloc[:, is_skewed], nrows=2, figsize=(10,5))
featureplot(dfnum_t2.iloc[:, ~is_skewed], nrows=2, figsize=(10, 5))
Explanation: Now our originally skewed features look more symmetric.
End of explanation
dfnum_t2.index.names = ['Dataset', 'Id']
dfnum_t2.head()
dfnum_t2.to_csv('transformed_dataset_dfnum_t2.csv', index=True)
nfeatures = dfnum_t2.columns
target_t = logtarget.apply(standardize)
target_t.head()
Explanation: Save transformed numerical data
Use the storage magic to communicate between notebooks.
End of explanation
dfnum_t2.head()
corr = pd.DataFrame(data=dfnum_t2.loc['train',:].apply(lambda feature: sp.stats.pearsonr(feature, target_t['SalePrice'])), 
                    columns=['pearsonr'])
corr = corr.assign(correlation=corr.applymap(lambda x: x[0]),
                pvalue=corr.applymap(lambda x: x[1]))
corr = corr.drop('pearsonr', axis=1)
corr.head()
corr.sort_values('pvalue', ascending=False)['correlation'].plot(kind='barh')
corr.sort_values('pvalue').head()
corr.sort_values('pvalue').tail()
Explanation: Feature selection
We're now in a good position to identify the key numerical features. Those should be hightly correlated with the sale price.
End of explanation
min_correlation = 0.2
key_features = corr[np.abs(corr['correlation'] > min_correlation)].sort_values(by='correlation', ascending=False).index.values
key_features, key_features.size
%store key_features
Explanation: Let's keep only the features that have a high enough correlation with the price (correlation less than 0.2)
End of explanation
<END_TASK>