<jupyter_start><jupyter_text>Imports<jupyter_code>%config InlineBackend.figure_format = 'retina'
%matplotlib inline

import numpy as np
import scipy as sp
from scipy import stats
import matplotlib.pyplot as plt
import statsmodels.api as sm
import pandas as pd
import seaborn as sns
import util

pd.options.display.max_columns=25
sns.set_style('white')<jupyter_output>/Users/scott/anaconda/lib/python3.6/site-packages/statsmodels/compat/pandas.py:56: FutureWarning: The pandas.core.datetools module is deprecated and will be removed in a future version. Please use the pandas.tseries module instead.
  from pandas.core import datetools<jupyter_text>Load data<jupyter_code>df, df_shops, df_ingredients = util.load_burritos()
df.tail(20)<jupyter_output><empty_output><jupyter_text>Analyze correlations between features<jupyter_code># Throw out outlier data
df_use = df.copy()
df_use = df_use[df_use['Cost'] < 9]
df_use = df_use[df_use['Cost'] > 5]
df_use = df_use[np.logical_or(df_use['Volume'] < 1.5, df_use['Volume'].isnull())]
df_use = df_use[np.logical_or(df_use['Volume'] > 0.5, df_use['Volume'].isnull())]

# Compute correlation matrix
m_corr = ['Hunger','Cost','Volume','Tortilla','Temp','Meat','Fillings','Meat:filling',
               'Uniformity','Salsa','Wrap','overall']
dfcorr = df_use[m_corr].corr(method='spearman')
for feat in m_corr:
    dfcorr.loc[feat, feat] = np.nan

# Determine the most highly correlated feature for each other feature
for feat in m_corr:
    feat2 = dfcorr[feat].argmax()
    print(feat, feat2,
          stats.spearmanr(df_use[feat], df_use[feat2], nan_policy='omit')[0],
          stats.spearmanr(df_use[feat], df_use[feat2], nan_policy='omit')[1])
# Compute correlation matrix for restaurants
df_rest = df_use.groupby('Location')[m_corr].mean()
df_rest_corr = df_rest.corr(method='spearman')
for feat in m_corr:
    df_rest_corr.loc[feat, feat] = np.nan

# Determine the most highly correlated feature for each other feature
for feat in m_corr:
    feat2 = df_rest_corr[feat].argmax()
    print(feat, feat2,
          stats.spearmanr(df_rest[feat], df_rest[feat2], nan_policy='omit')[0],
          stats.spearmanr(df_rest[feat], df_rest[feat2], nan_policy='omit')[1])
print(stats.spearmanr(df_rest['Meat'], df_rest['Meat:filling'], nan_policy='omit'))
print(stats.spearmanr(df_use['Meat'], df_use['Meat:filling'], nan_policy='omit'))

print(stats.spearmanr(df_rest['Fillings'], df_rest['Meat:filling'], nan_policy='omit'))
print(stats.spearmanr(df_use['Fillings'], df_use['Meat:filling'], nan_policy='omit'))

print(stats.spearmanr(df_rest['Fillings'], df_rest['Meat'], nan_policy='omit'))
print(stats.spearmanr(df_use['Fillings'], df_use['Meat'], nan_policy='omit'))
# within taco shop studies
most_places = df_use['Location'].value_counts().index[:5]
for r in most_places:
    df_temp = df_use[df_use['Location']==r]
    print(r, stats.spearmanr(df_temp['Fillings'], df_temp['Meat'], nan_policy='omit'))
plt.figure(figsize=(5,5))
plt.plot(df_rest['Meat'], df_rest['overall'], 'k.')
plt.xlim((1,5.1))
plt.ylim((1,5.1))
plt.plot([1,5], [1,5], 'k--')
plt.xlabel('Average meat rating', size=15)
plt.ylabel('Average overall rating', size=15)

plt.figure(figsize=(5,5))
plt.plot(df_rest['Temp'], df_rest['overall'], 'k.')
plt.xlim((.9,5.1))
plt.ylim((.9,5.1))
plt.plot([1,5], [1,5], 'k--')
plt.xlabel('Average temperature rating', size=15)
plt.ylabel('Average overall rating', size=15)<jupyter_output><empty_output><jupyter_text>dendrogram to cluster burritos<jupyter_code>plt.figure(figsize=(10,20))
m_corr = ['Hunger','Tortilla','Temp','Meat','Fillings','Meat:filling','Uniformity','Salsa','Wrap','overall']
sns.clustermap(df_use[m_corr].dropna())
plt.figure(figsize=(10,20))
g = sns.clustermap(df_rest.corr(method='spearman'), vmin=0, vmax=1, cmap='mako')<jupyter_output><empty_output><jupyter_text>Hierarchy scipy dendrogram<jupyter_code>import scipy.cluster.hierarchy as hier

Y = hier.linkage(df_rest.corr(method='spearman'),
                 method='complete',
                 metric='cosine')

plt.figure(figsize=(3,15))

Z = hier.dendrogram(Y, orientation='left',
                    labels=df_rest.columns,
                    leaf_font_size=12)<jupyter_output><empty_output><jupyter_text>My ratings<jupyter_code># Compare GLMs for me and others
m_lm = ['Hunger','Tortilla','Temp','Meat','Fillings','Meat:filling',
        'Uniformity','Salsa','Wrap','Cost']

df_scott = df_use[df_use['Reviewer']=='scott']
df_emily = df_use[df_use['Reviewer']=='emily']
df_ricardo = df_use[df_use['Reviewer']=='ricardo']
df_notscott = df_use[df_use['Reviewer']!='scott']

# Remove incomplete data 
dffull = df_scott[np.hstack((m_lm,'overall'))].dropna()

# Run model and display results
X = sm.add_constant(dffull[m_lm])
y = dffull['overall']
my_glm = sm.GLM(y,X)
res = my_glm.fit()
print('Explained variance: {:.1f}%'.format(100*(1 - np.var(res.resid_pearson) / np.var(y))))
print(res.summary())

# Remove incomplete data 
dffull = df_notscott[np.hstack((m_lm,'overall'))].dropna()

# Run model and display results
X = sm.add_constant(dffull[m_lm])
y = dffull['overall']
my_glm = sm.GLM(y,X)
res = my_glm.fit()
print('Explained variance: {:.1f}%'.format(100*(1 - np.var(res.resid_pearson) / np.var(y))))
print(res.summary())

# Remove incomplete data 
dffull = df_emily[np.hstack((m_lm,'overall'))].dropna()

# Run model and display results
X = sm.add_constant(dffull[m_lm])
y = dffull['overall']
my_glm = sm.GLM(y,X)
res = my_glm.fit()
print('Explained variance: {:.1f}%'.format(100*(1 - np.var(res.resid_pearson) / np.var(y))))
print(res.summary())

# Remove incomplete data 
dffull = df_ricardo[np.hstack((m_lm,'overall'))].dropna()

# Run model and display results
X = sm.add_constant(dffull[m_lm])
y = dffull['overall']
my_glm = sm.GLM(y,X)
res = my_glm.fit()
print('Explained variance: {:.1f}%'.format(100*(1 - np.var(res.resid_pearson) / np.var(y))))
print(res.summary())
df_use['Reviewer'].value_counts()
df_use['Location'].value_counts()<jupyter_output><empty_output><jupyter_text>Compare restaurants in each dimension<jupyter_code>sns.set_style('whitegrid')
restaurants_compare = ['taco stand', 'california burritos', "rigoberto's taco shop", 'taco villa', 'lucha libre north park']

m_lm = ['Tortilla','Temp','Meat','Fillings','Meat:filling','Uniformity','Salsa','Synergy','Wrap','overall']
df_plt = df_use[df_use['Location'].isin(restaurants_compare)]
df_plt = df_plt[np.hstack((m_lm,'Location'))]
df_plt.set_index('Location', inplace=True)
df_plt = df_plt.stack().reset_index().rename(columns={'level_1': 'feature', 0: 'value'})

plt.figure(figsize=(16,7))
sns.barplot(x='feature', y='value', hue='Location', data=df_plt, hue_order=restaurants_compare)
plt.gca().legend(title='Location', fontsize=15, loc='center left', bbox_to_anchor=(1, 0.5))
plt.setp(plt.gca().get_legend().get_title(), fontsize='15')
plt.xticks(size=15)
plt.yticks(np.arange(2.5,5,.5), size=15)
plt.xlabel('Burrito feature', size=15)
plt.ylabel('Average value (out of 5)', size=15)
plt.ylim((2.5,4.6))


m_lm = ['Cost']
df_plt = df_use[df_use['Location'].isin(restaurants_compare)]
df_plt['dummy'] = 1

plt.figure(figsize=(3,7))
sns.barplot(x='dummy', y='Cost', hue='Location', data=df_plt, hue_order=restaurants_compare)
plt.gca().legend(title='Location', fontsize=15, loc='center left', bbox_to_anchor=(1, 0.5))
plt.setp(plt.gca().get_legend().get_title(), fontsize='15')
plt.xticks([])
plt.yticks(np.arange(0,9), size=15)
plt.xlabel('')
plt.ylabel('Average cost ($)', size=15)
plt.ylim((0,8))

m_lm = ['Volume']
df_plt = df_use[df_use['Location'].isin(restaurants_compare)]
df_plt['dummy'] = 1

plt.figure(figsize=(3,7))
sns.barplot(x='dummy', y='Volume', hue='Location', data=df_plt, hue_order=restaurants_compare)
plt.gca().legend(title='Location', fontsize=15, loc='center left', bbox_to_anchor=(1, 0.5))
plt.setp(plt.gca().get_legend().get_title(), fontsize='15')
plt.xticks([])
plt.yticks(np.arange(0,1.1,.2), size=15)
plt.xlabel('')
plt.ylabel('Average size (L)', size=15)
plt.ylim((0, 1))<jupyter_output>/Users/scott/anaconda/lib/python3.6/site-packages/ipykernel_launcher.py:23: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy
/Users/scott/anaconda/lib/python3.6/site-packages/ipykernel_launcher.py:37: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy<jupyter_text>Most distinctive featureDifference between across-restaurant variance and within-restaurant variance<jupyter_code>feats = ['Tortilla','Temp','Meat','Fillings','Meat:filling','Uniformity','Salsa','Synergy','Wrap','overall', 'Hunger']

# Compute variance within restaurants
var_within_rest = df_use.groupby('Location')[feats].var().mean()
var_across_rest = df_use.groupby('Location')[feats].mean().var()
var_within_by_across = var_within_rest / var_across_rest

df_var = pd.DataFrame({'within shop': var_within_rest,
                       'across shop': var_across_rest, 
                       'within : across': var_within_by_across})
df_var = df_var.stack().reset_index().rename(columns={'level_0':'Feature', 'level_1': 'type', 0: 'Variance'})

plt.figure(figsize=(10,5))
sns.barplot(x='type', y='Variance', hue='Feature', data=df_var, palette='viridis',
            hue_order=var_within_by_across.sort_values().index,
            order = ['within shop', 'across shop', 'within : across'])
plt.gca().legend(title='Location', fontsize=15, loc='center left', bbox_to_anchor=(1, 0.5))
plt.setp(plt.gca().get_legend().get_title(), fontsize='15')
plt.xlabel('Variance type', size=15)
plt.ylabel('Normalized Variance', size=15)
plt.xticks(size=15)
plt.yticks(np.arange(0, 2.5, .5), size=15)
plt.ylim((0,2))

feats = ['Volume', 'Cost']

# Compute variance within restaurants
var_within_rest = df_use.groupby('Location')[feats].var().mean()
var_across_rest = df_use.groupby('Location')[feats].mean().var()
var_within_by_across = var_within_rest / var_across_rest
df_var = pd.DataFrame({'within shop': var_within_rest,
                       'across shop': var_across_rest, 
                       'within : across': var_within_by_across})
df_var.head()<jupyter_output><empty_output><jupyter_text>Distribution of features<jupyter_code>cols_plot = ['Tortilla', 'Temp', 'Meat', 'Fillings', 'Meat:filling', 'Uniformity', 'Salsa', 'Synergy', 'Wrap']
df_plt = df[cols_plot]
df_plt = df_plt.stack().reset_index().drop('level_0', axis=1).rename(columns={'level_1': 'feature', 0: 'value'})

plt.figure(figsize=(6,6))
sns.violinplot(y='feature', x='value', data=df_plt)
plt.xlim((1,5))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18215069495245062, 'fraction_numerical': 0.018653986832479885, 'mean_word_length': 6.718419195483416, 'pattern_counts': {'":': 0, '<': 37, '<?xml version=': 0, '>': 37, 'https://': 0, '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': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '644163', 'license': 'MIT', 'n_tokens_mistral': 4440, 'n_tokens_neox': 4127, 'n_words': 674, 'path': 'April2018.ipynb', 'repo_name': 'srcole/burritos', 'size': '10936'}