<jupyter_start><jupyter_text>Table of Content [1. Reading the data](0)[2. data manipulation](1)&nbsp;&nbsp;&nbsp;&nbsp;[2.1 first contact with the data](2)&nbsp;&nbsp;&nbsp;&nbsp;[2.2 accessing specific parts of the data](3)&nbsp;&nbsp;&nbsp;&nbsp;[2.3 adding and combining columns](4)&nbsp;&nbsp;&nbsp;&nbsp;[Exercises 01 - data manipulation](5)[3. data description and representation](6)&nbsp;&nbsp;&nbsp;&nbsp;[3.1 basic description - common summary statistics](7)&nbsp;&nbsp;&nbsp;&nbsp;[3.2 Representing one column - histograms and density line](8)&nbsp;&nbsp;&nbsp;&nbsp;[Exercises 02 - histograms](9)&nbsp;&nbsp;&nbsp;&nbsp;[3.3 accounting for categories in the data](10)&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;[3.3.1 Boxplots](11)&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;[3.3.2 violin-plots](12)&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;[3.3.3 barplots and their error bars](13)&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;[3.3.4 showing several categorizations at the same time](14)&nbsp;&nbsp;&nbsp;&nbsp;[Exercises 03 - representing categories](15)[4 Additional types of plots and some *caveat* about graphical representations](16)  importing, manipulating, and representing data The basis of any statistical analysis is the underlying data.A data-set is typically presented as a file containing information formatted as a table: * each line correspond to an observation ( individual, sample, ... ) * each column correspond to a measured variable ( height, sex, gene expression, ... )To read data file and manipulate the data, we will rely on [pandas](https://pandas.pydata.org/)Pandas is a "high-level" module, designed for statistics/exploratory analysis.A great strength of pandas is its **DataFrame** which emulates many of the convenient behavior and syntax of their eponym counterpart in the **R** language.To graphically represent the data, we will rely on [seaborn](https://seaborn.pydata.org/index.html).Seaborn is designed to work hand-in-hand with pandas DataDrame to produce **efficient data representation** from fairly simple commands. They propose very good tutorials as well as a gallery (with associated code) that can get you started quickly.<jupyter_code>import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
print('OK')<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) 1. Reading the data what is the file name? location?what is the saprator between fields??`pd.read_table` is a generalistic function to read table. Aside from the name of the file to read, here are some useful parameters:* `sep` : separator between columns (by default '\t')* `header` : Row number(s) to use as the column names. By default it will use the first line as a header. use `header=None` if the file does not contain column name.* `skiprows` : Line numbers to skip (0-indexed) or number of lines to skip (int) at the start of the file.* (`true_values`/`false_values` can also be useful)Of course you can learn (much) more using `help(pd.read_table)`.The data we will be playing with is census data of the swiss population for the year 1880 [source](https://www.bfs.admin.ch/bfs/en/home/statistics/population/surveys/census.assetdetail.11047857.html).The data-file is named `"swiss_census_1880.csv"` and like its extension suggests, it contains **C**omma-**S**eparated **V**alues.<jupyter_code>import pandas as pd

df = pd.read_table("data/swiss_census_1880.csv" , sep=',') 
#try to see what happens when sep has a different value

df.head() # this returns the 5 first lines of the table<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) 2. data manipulation [back to the toc](toc) 2.1 first contact with the data Gathering basic information about the data-set is fairly easy :<jupyter_code>numberRows , numberCols = df.shape
print('rows:',numberRows , 'columns:',numberCols)
print('column names')
df.columns
print( "columns types:\n" ,df.dtypes )# lists the types<jupyter_output><empty_output><jupyter_text>Here we can see that most values are integers, save for names (town names, canton, ...). [back to the toc](toc) 2.2 accessing specific parts of the data One can access a column just by using `df[]` :<jupyter_code>df['Foreigner']<jupyter_output><empty_output><jupyter_text>To access a **specific section of the data, use the `.loc[ : , : ]` operator.**It takes two slicing arguments thatt, respectively, targets rows and columns. A `:` selects all elements.<jupyter_code>df.loc[ 0:3 , : ] ### first 4 rows, all columns
df.loc[ 0:3 , ['town name' , 'Total'] ] ### first 4 rows, columns 'town name' and 'Total'
df.loc[ : , 'Total' ] ### all rows, column 'Total'<jupyter_output><empty_output><jupyter_text>> There also exists the `.iloc[ : , : ]` operator to select columns (and rows) by index instead of names. See this [pandas tutorial](https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html) for more details and recipes. Another, powerful, way of accessing specific part of the data is by defining a mask, which will **filter the data through a particular condition**.<jupyter_code>maskVD = df.loc[ : , 'canton' ] == 'VD'  
maskVD.value_counts()<jupyter_output><empty_output><jupyter_text>The mask (in fact a `pandas.core.series.Series`), is in effect a list of values that are ```True``` or `False` depending on whether or not they satisfy the defined condition (canton is equal to 'VD', here).The mask can then be applied to the `DataFrame`.<jupyter_code>df.loc[ maskVD  , ['town name','canton','canton name'] ]<jupyter_output><empty_output><jupyter_text>Masks may be combined to produce more complex selection criteria, using **`&`** (logical and) and **`|`** (logical or).<jupyter_code># towns in the canton of Zurich and with more than 10 000 registered inhabitant
mask = (df.loc[ : , 'canton' ] == 'ZH') & ( df.loc[ : , 'Total' ] >= 10000  ) 

df.loc[ mask , ]<jupyter_output><empty_output><jupyter_text>**micro-exercise :** Select towns with less than 1000 inhabitants (column `Total`), (*optional*: display only town name and number of inhabitants) [back to the toc](toc) 2.3 adding and combining columns New columns are easy to create from existing ones :<jupyter_code>df['14+ y.o.'] = df['15-59 y.o.'] + df['60+ y.o.']
df['14+ y.o.'].head()<jupyter_output><empty_output><jupyter_text>Removing columns is about as easy:<jupyter_code>df = df.drop(columns='14+ y.o.') # use the 'rows' argument to remove rows instead
print("is '14+ y.o.' part of the columns : " , '14+ y.o.' in df.columns)<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) Exercises 01 - data manipulation 1. Select towns with less than 1000 inhabitants, or with more than 1 Foreigner. How many such towns are there? (*hint:* build from the last micro-exercise)2. Create a new column is the `DataFrame` representing the fraction of population which is Reformed in each town.(*optional : What is the minimum/maximum value for this fraction?*)<jupyter_code># %load solutions/solution_01_01.py<jupyter_output><empty_output><jupyter_text>InterludeFor the rest of the notebook, we will be not be using the counts themselves much, but rather their corresponding fraction of the registered populations. Le'ts create another DataFrame to ease our process:<jupyter_code>dfFractions = df.loc[:, ['town name','canton','canton name' , 'Total'] ]
columnsOfInterest = [ 'Swiss', 'Foreigner','Male', 'Female', '0-14 y.o.', '15-59 y.o.', '60+ y.o.', 'Reformed','Catholic', 'Other', 'German speakers', 'Franch speakers','Italian speakers', 'Romansche speakers','Non-national tongue speakers' ]
for x in columnsOfInterest:
    dfFractions[x] = df[x]/df['Total']
dfFractions.head()<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) 3. data description and representation  [back to the toc](toc) 3.1 basic description - common summary statistics<jupyter_code>dfFractions.describe()<jupyter_output><empty_output><jupyter_text>The `.describe()` function gives information about all numerical columns in the data-set at once (note that non-numerical values are absent here).It is very useful not only to get a first impression on the data-set, but also to catch eventual errors in the data : a negative number where there should be only positive values, missing values (NAs), ...**Should I add some problematic values to the data in order to demonstrate?** `.describe()` gives access to some of the most commonly used summary statistics:* (arithmetic) **mean** : ${\displaystyle \bar{x}={\frac {1}{n}}\sum _{i=1}^{n}x_{i}}$  or, for coders : `sum(x)/len(x)`* **standard deviation** (a.k.a. std, stdev) : this corresponds to the average of the absolute difference to the mean. It is the **square root of the variance**.* **minimum** and **maximum**: smallest and biggest value among the data. Looking at them can help detect outliers.* **Quartiles** : they correspond to the value such that    * 25% (first quartile, Q1),     * 50% (second quartile, Q2, median), or    * 75% (second quartile, Q3)          of the values are lower than them. They are less sensitive than the mean to outlier values.To these, one could also add the **mode**: the most common value (discrete data), or the value with the maximum density (continuous data). Data distributions with multiple peaks (e.g., in the presence of multiple subpopulations) are called multimodal.Sometimes mode, median and mean coincide, other times they can be quite different. This is best demonstrated demonstrated visually, which is why we will be talking about ... [back to the toc](toc) 3.2 Representing one column - histograms and density line The simplest way to represent the distribution of a set of values using `seaborn` is to use the function `sns.histplot`:<jupyter_code>sns.histplot(dfFractions['Swiss'] )<jupyter_output><empty_output><jupyter_text>You can see an histogram (light-blue).The *x-axis* corresponds to the fraction of Swiss people, and the *y-axis* corresponds the **count**, which for the histogram corresponds to the number of values falling in a given *bins* (i.e. a bar of the histogram).As with any graphical function, you can customize almost everything.Let's add a density line and change the color :<jupyter_code>sns.histplot(dfFractions['Swiss'], kde=True , color='red' ).set_title('histogram alone')<jupyter_output><empty_output><jupyter_text>> BTW, that is how you set-up colors and add a title on seaborn**Adding the density line switch the y-axis from count to density.**A different number of *bins* shows different aspects of the distribution of your data and so it is important to choose one that is appropriate to your sample size and data. By default, seaborn automatically infers the number of *bins*. You may stick with the original or change it :<jupyter_code># Set up a figure with multiple panels, here a 2 by 2 grid
f, axes = plt.subplots( 2, 2, figsize=(7, 7) )

# then you can specify where each plot goes on the figure with the ax argument ([0,0] is the top left corner)

# Plot a simple histogram with binsize determined automatically
sns.histplot(dfFractions['0-14 y.o.'], kde=False, color="b", ax=axes[0, 0]).set_title('automatic')

# Plot a simple histogram with binsize 5, 10 , 1000
sns.histplot(dfFractions['0-14 y.o.'], bins=5   , kde=False, color="b", ax=axes[0, 1]).set_title('5 bins')
sns.histplot(dfFractions['0-14 y.o.'], bins=10  , kde=False, color="b", ax=axes[1, 0]).set_title('10 bins')
sns.histplot(dfFractions['0-14 y.o.'], bins=1000 , kde=False, color="b", ax=axes[1, 1]).set_title('1000 bins')

plt.tight_layout()# this makes the panels margins and setup more graceful
sns.kdeplot(dfFractions['Reformed'], color="y")<jupyter_output><empty_output><jupyter_text>Now that we have seen how to represent the distribution of values in a column, let's come back to our summary statistics, in particular the **mean**, the **median** and the **mode**<jupyter_code>import seaborn as sns
import matplotlib.pyplot as plt

# Here we just define a small function for plotting a distribution with the mean median and mode 
def plotWithMeanMedianMode( dat , ax):

    mode=dat.mode()[0] #we only select the first mode
    mean=dat.mean()
    median=dat.median()

    sns.histplot( dat , ax=ax , kde=True) # line for histogram and density line

    ax.axvline(mean, color='r', linestyle='--' )
    ax.axvline(median, color='g', linestyle='-')
    ax.axvline(mode, color='b', linestyle='-')
    ax.legend({'Mean':mean,'Median':median,'Mode':mode})
f, axes = plt.subplots( 3, 1, figsize=(15, 10) )

plotWithMeanMedianMode( dfFractions['0-14 y.o.'] , ax=axes[0])
plotWithMeanMedianMode( dfFractions['Foreigner'] , ax=axes[1])
plotWithMeanMedianMode( dfFractions['Reformed'] , ax=axes[2])
plt.show()<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) Exercises 02 - histograms 1. plot the distribution of the total number of habitants. Try to choose an appropriate mode of representation (histogram, density line? number of bins?)2. try to call `sns.histplot` twice in a row, once with to plot the fraction of Foreigner and the other for the fraction of Swiss. What happens?3. plot the distribution of the fraction of catholics in the canton of Zurich.<jupyter_code># %load -r 1-7 solutions/solution_01_02.py
# %load -r 10-14 solutions/solution_01_02.py
# %load -r 15- solutions/solution_01_02.py<jupyter_output><empty_output><jupyter_text>Of course, we have only scratched the surface of this form of representation. You can learn more from the [seaborn official tutorial](https://seaborn.pydata.org/tutorial/distributions.htmlplotting-univariate-distributions) [back to the toc](toc) 3.3 accounting for categories in the data  Being able to show the distribution of the fraction of people aged 14 years old or less for your whole data-set(here, at national level) is all good, but often we want to compute and represent this information for different categories.Here the categories would be the canton, but depending on data-sets they could be a genotype (Wild-Type versus KO), sex, experimental condition, ...Of course, we could manually compute the means and plot an histogram for each category in our data-set, but pandas and seaborn offer much more efficient routines for this. We will also create a column that describes the main religion and main languague for each town:<jupyter_code>dfFractions['majority religion'] = dfFractions.loc[ : , ['Reformed','Catholic','Other'] ].idxmax(axis=1)
dfFractions['majority language'] = dfFractions.loc[ : , ['German speakers', 
                                                         'Franch speakers',
                                                         'Italian speakers', 
                                                         'Romansche speakers',
                                                         'Non-national tongue speakers'] ].idxmax(axis=1)<jupyter_output><empty_output><jupyter_text>`.idxmax(axis=1)` picks out the column with maximum value for each row. Computing a summary statistics on a pandas DataFrame is done **using the `df.groupby(...)` method and the applying some function** to the result of that grouping.<jupyter_code>grouped = df.groupby('canton name')## grouping by canton
grouped['Total'].sum() ## total number of registered inhabitant for each canton 
dfFractions.groupby('majority language')['Catholic'].mean() ## mean fraction of caholics in towns depending on the majority language
def getSmallestTown( data ):
    ''' returns the name and population of the town with minimal number of inhabitants '''
    indexsmallestTown = data['Total'].idxmin()
    return data['town name'][indexsmallestTown] , data['Total'][indexsmallestTown]

grouped.apply(getSmallestTown) ## name and population of the town with minimal number of inhabitants for each canton<jupyter_output><empty_output><jupyter_text>Computing these summary statistics is all good, but you can also vizualize them using **seaborn master function sns.catplot**<jupyter_code>sns.catplot( x = 'majority language' , y='Catholic' ,  data=dfFractions)<jupyter_output><empty_output><jupyter_text>The `kind` argument let's you control the overall look of the plot. IIt can be:* 'strip' : this is the default, showing all data points. * **'box'** : the famous boxplot* **'violin'** : an alternative to the boxplot using density lines* **'bar'** : the (in)famous barplot* 'swarm' : similar to 'strip' but with another way of aranging the points* 'boxen' : some intermediary between a boxplot and a violin plot* 'point' : alternative to the barplot where only the top point is shown<jupyter_code>kinds = ['box','violin','bar','boxen','strip','point']#,'swarm'] # swarm takes a long time to compute


for i,k in enumerate(kinds):
    sns.catplot( x = 'Male' , y='majority religion' ,  data=dfFractions , kind = k ,height=2, aspect=5 )

# note 2 arguments here :
# * height : height of the plot
# * aspect : widht/height ratio. -> higher aspect = wider figure<jupyter_output><empty_output><jupyter_text>We are going to concentrate on **boxplot**, **violin-plot** and **barplot** as they arguably the more used (and misused) plot encountered in publications. [back to the toc](toc) 3.3.1 Boxplots Boxplots are made of :* a box : the first, second and third quartile of the data* whiskers : outliers But a picture is worth a thousand words:![boxplot](images/boxplot.png) * Q1/3 : first/third quartile (25%/75% of the data below) * IQR : interquartile range : distance between the first and third quartile While they are very powerful, **boxplot have trouble with multimodal data**<jupyter_code>sns.catplot( x = 'German speakers' , y='majority religion' ,  data=dfFractions , kind = 'box' ,height=2, aspect=5 )
sns.catplot( x = 'German speakers' , y='majority religion' ,  data=dfFractions , kind = 'strip' ,height=2, aspect=5 )<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) 3.3.2 violin-plots These are representations takes boxplots one step further by drawing density lines instead of box and whiskers.<jupyter_code>sns.catplot( x = 'German speakers' , y='majority religion' ,  data=dfFractions , kind = 'violin' ,height=2, aspect=5 )<jupyter_output><empty_output><jupyter_text>They are generally more versatile than boxplots. However, the estimation of the density line means that they work better when the amount of data is large, and they give **strange results when the represented values are bounded** (see here, the violin plot show a density for fraction of german speakers under 0 or above 1.0 ...).You can solve this with the option `cut=0` (IMO they should make that one the default value).<jupyter_code>sns.catplot( x = 'German speakers' , y='majority religion' ,  data=dfFractions , kind = 'violin' ,height=2, aspect=5 , cut=0)<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) 3.3.3 barplots and their error bars This representation is widely used, but it is, demonstrably, a less than optimal way of representing a distribution. At its code, it is a representation of two summary statistics:* the **mean**, which is the size of the bar* an **error value** which can be:    * range : bars extend from smallest to largest value    * standard deviation ($sigma$): square root of the variance    * standard error of the mean (SEM): $sem = \frac{\sigma}{\sqrt{n}}$    * 95% confidence interval (CI): $1.96 * \frac{\sigma}{\sqrt{n}}$The standard deviation is a measure of the variability in the data, while the **standard error of the mean is a measure of the precision of the estimate of the mean**.<jupyter_code>sns.catplot( x = 'German speakers' , y='majority religion' , data=dfFractions , height=2, aspect=5,
            kind = 'bar' , ci='sd').set(title='standard deviation')
sns.catplot( x = 'German speakers' , y='majority religion' , data=dfFractions , height=2, aspect=5,
            kind = 'bar' , ci=95 ).set(title='95% confidence interval')<jupyter_output><empty_output><jupyter_text>**You should ALWAYS specify which type of error bar you are using.** Given all this, one can reasonnably ask the question : **Is there a *best* way ?**While the answer is not clear about the best way to represent a distribution over multiple categories, we definitely have some good answers regarding bad/misleading ways : ![barplotnonsense](images/barplotNonsense.png) From https://stekhoven.shinyapps.io/barplotNonsense/Barplots :* Contain the same amount of information as a purely numeric summary (mean ± error)* **Avoid it if possible*** Always indicate what type of error bars you are using* Many articles have been written on the subject of bar charts and possible alternatives, e.g. [a Nature Biomedical Engineering Editorial](https://www.nature.com/articles/s41551-017-0079), [<NAME> (2014)](https://www.nature.com/articles/nmeth.2807), [a Nature Methods Editorial](https://www.nature.com/articles/nmeth.2837), [<NAME> (2011)](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3087125/) and [Weissgerber et al (2015)](https://journals.plos.org/plosbiology/article?id=10.1371/journal.pbio.1002128) [back to the toc](toc) 3.3.4 showing several categorizations at the same time `seaborn` offers an easy way to add categories to an existing plot with the `hue` argument:<jupyter_code>#example with a density line plot

k = 'box' 

category='majority religion'
represented_variable="Italian speakers" 

sns.kdeplot( x = represented_variable , #hue= category , # uncomment previous to see the effect of hue
             data=dfFractions )
kinds = ['box','violin','bar','boxen','strip','point']#,'swarm'] # swarm takes a long time to compute

categorizations=["majority language","majority religion"]

## uncomment this to switch the categories
#categorizations = categorizations[::-1]

represented_variable="Total" 

for i,k in enumerate(kinds):
    g=sns.catplot( x = represented_variable , y= categorizations[0] , hue= categorizations[1] ,  
                data=dfFractions , kind = k , orient='h',height=2, aspect=5 )
    g.set(xscale="log")<jupyter_output><empty_output><jupyter_text>[back to the toc](toc) Exercises 03 - representing categories Represent the proportion of people more 60 years old (`'60+ y.o.'`) across all cantons Choose the most appropriate kind of plot.<jupyter_code># %load  solutions/solution_01_03.py<jupyter_output><empty_output><jupyter_text>interlude - writing data and plotTo write your data-set to a file, use :<jupyter_code>dfFractions.to_csv('census1880_fractions.csv')<jupyter_output><empty_output><jupyter_text>and for figures :<jupyter_code># make you plot as usual, redirect them to a variable
represented_variable , category = 'Total' , 'majority language'
myPlot = sns.catplot( x = represented_variable , y= category ,
             data=dfFractions , kind = k , orient='h',height=10, aspect=2 )

myPlot.savefig("output.png")
# many formats can be specified with the name of the file
# for example : .pdf, .svg, .jpg, ...
# for raster formats, the dpi argument can help define the image quality
# For multipanel figures :
f, axes = plt.subplots( 2, 2, figsize=(7, 7) )

sns.histplot(dfFractions['0-14 y.o.'], kde=False, color="b", ax=axes[0, 0]).set_title('automatic')
sns.histplot(dfFractions['0-14 y.o.'], bins=5   , kde=False, color="b", ax=axes[0, 1]).set_title('5 bins')
sns.histplot(dfFractions['0-14 y.o.'], bins=10  , kde=False, color="b", ax=axes[1, 0]).set_title('10 bins')
sns.histplot(dfFractions['0-14 y.o.'], bins=1000 , kde=False, color="b", ax=axes[1, 1]).set_title('1000 bins')
plt.tight_layout()# this makes the panels margins and setup more graceful

f.savefig('output_multipanel.pdf')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11840095060261416, 'fraction_numerical': 0.017611610931930064, 'mean_word_length': 5.13671875, 'pattern_counts': {'":': 0, '<': 123, '<?xml version=': 0, '>': 127, 'https://': 11, 'lorem ipsum': 0, 'www.': 5, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '30', 'dir': 'jupyter-structured-clean-dedup', 'id': '262725', 'license': 'CC0-1.0', 'n_tokens_mistral': 7464, 'n_tokens_neox': 6909, 'n_words': 2663, 'path': 'Intro_to_stats_in_life_sciences/01_data_manipulation_and_representation.ipynb', 'repo_name': 'SebastienBoyer/Teaching', 'size': '23587'}