π day special: How to use Bioconductor to find empirical evidence in support of π being a normal number

Editor’s note: Today 3/14/15 at some point between  9:26:53 and 9:26:54 it was the most π day of them all. Below is a repost from last year. 

Happy π day everybody!

I wanted to write some simple code (included below) to the test parallelization capabilities of my  new cluster. So, in honor of  π day, I decided to check for evidence that π is a normal number. A normal number is a real number whose infinite sequence of digits has the property that picking any given random m digit pattern is 10^−m. For example, using the Poisson approximation, we can predict that the pattern “123456789” should show up between 0 and 3 times in the first billion digits of π (it actually shows up twice starting, at the 523,551,502-th and  773,349,079-th decimal places).

To test our hypothesis, let Y₁, …, Y₁₀₀ be the number of “00”, “01”, …,”99” in the first billion digits of  π. If  π is in fact normal then the Ys should be approximately IID binomials with N=1 billon and p=0.01.  In the qq-plot below I show Z-scores (Y - 10,000,000) /  √9,900,000) which appear to follow a normal distribution as predicted by our hypothesis. Further evidence for π being normal is provided by repeating this experiment for 3,4,5,6, and 7 digit patterns (for 5,6 and 7 I sampled 10,000 patterns). Note that we can perform a chi-square test for the uniform distribution as well. For patterns of size 1,2,3 the p-values were 0.84, 0.89, 0.92, and 0.99.

Another test we can perform is to divide the 1 billion digits into 100,000 non-overlapping segments of length 10,000. The vector of counts for any given pattern should also be binomial. Below I also include these qq-plots.

These observed counts should also be independent, and to explore this we can look at autocorrelation plots:

To do this in about an hour and with just a few lines of code (included below), I used the Bioconductor Biostrings package to match strings and the foreach function to parallelize.

`Editor’s note: Today 3/14/15 at some point between  9:26:53 and 9:26:54 it was the most π day of them all. Below is a repost from last year. 

Happy π day everybody!

I wanted to write some simple code (included below) to the test parallelization capabilities of my  new cluster. So, in honor of  π day, I decided to check for evidence that π is a normal number. A normal number is a real number whose infinite sequence of digits has the property that picking any given random m digit pattern is 10^−m. For example, using the Poisson approximation, we can predict that the pattern “123456789” should show up between 0 and 3 times in the first billion digits of π (it actually shows up twice starting, at the 523,551,502-th and  773,349,079-th decimal places).

To test our hypothesis, let Y₁, …, Y₁₀₀ be the number of “00”, “01”, …,”99” in the first billion digits of  π. If  π is in fact normal then the Ys should be approximately IID binomials with N=1 billon and p=0.01.  In the qq-plot below I show Z-scores (Y - 10,000,000) /  √9,900,000) which appear to follow a normal distribution as predicted by our hypothesis. Further evidence for π being normal is provided by repeating this experiment for 3,4,5,6, and 7 digit patterns (for 5,6 and 7 I sampled 10,000 patterns). Note that we can perform a chi-square test for the uniform distribution as well. For patterns of size 1,2,3 the p-values were 0.84, 0.89, 0.92, and 0.99.

Another test we can perform is to divide the 1 billion digits into 100,000 non-overlapping segments of length 10,000. The vector of counts for any given pattern should also be binomial. Below I also include these qq-plots.

These observed counts should also be independent, and to explore this we can look at autocorrelation plots:

To do this in about an hour and with just a few lines of code (included below), I used the Bioconductor Biostrings package to match strings and the foreach function to parallelize.

`

NB: A normal number has the above stated property in any base. The examples above a for base 10.

De-weaponizing reproducibility

1. The tools for reproducibility already exist, the barrier isn’t tools
2. We need to de-weaponize reproducibility
3. Prevention is the right approach to reproducibility

In terms of the first point, tools like iPython, knitr, and Galaxy can be used to all but the absolutely largest analysis reproducible right now.  Our group does this all the time with our papers and so do many others. The problem isn’t a lack of tools.

Speaking to point two, I think many people would agree that part of the issue is culture change. One issue that is increasingly concerning to me is the “weaponization” of reproducibility.  I have been noticing is that some of us (like me, my students, other folks at JHU, and lots of particularly junior computational people elsewhere) are trying really hard to be reproducible. Most of the time this results in really positive reactions from the community. But when a co-author of mine and I wrote that paper about the science-wise false discovery rate, one of the discussants used our code (great), improved on it (great), identified a bug (great), and then did his level best to humiliate us both in front of the editor and the general public because of that bug (not so great).

1.  Impose regulatory hurdles in the short term while people transition to reproducibility. One clear example of this is the Secret Science Reform Act which is a bill that imposes strict reproducibility conditions on all science before it can be used as evidence for regulation.
2. Humiliate people who aren’t good coders or who make mistakes in their code. This is what happened in my paper when I produced reproducible code for my analysis, but has also happened to other people.
3. Take advantage of people’s code to plagiarize/straight up steal work. I have stories about this I’d rather not put on the internet

Of the three, I feel like (1) and (2) are the most common. Plagiarism and scooping by theft I think are actually relatively rare based on my own anecdotal experience. But I think that the “weaponization” of reproducibility to block regulation or to humiliate folks who are new to computational sciences is more common than I’d like it to be. Until reproducibility is the standard for everyone - which I think is possible now and will happen as the culture changes -  the people who are the early adopters are at risk of being bludgeoned with their own reproducibility. As a community, if we want widespread reproducibility adoption we have to be ferocious about not allowing this to happen.

The elements of data analytic style - so much for a soft launch

Editor’s note: I wrote a book called Elements of Data Analytic Style. Buy it on Leanpub or Amazon! If you buy it on Leanpub, you get all updates (there are likely to be some) for free and you can pay what you want (including zero) but the author would be appreciative if you’d throw a little scratch his way. 

So uh, I was going to soft launch my new book The Elements of Data Analytic Style yesterday. I figured I’d just quietly email my Coursera courses to let them know I created a new reference. It turns out that that wasn’t very quiet. First this happened:

@jtleek @albertocairo @simplystats Instabuy. And apparently not just for me: it looks like you just Slashdotted @leanpub's website.

— Andrew Janke (@AndrewJanke) March 2, 2015

and sure enough the website was down:

then overnight it did something like 6,000+ units:

So lesson learned, there is no soft open with Coursera. Here is the post I was going to write though:

### Post I was gonna write

I have been doing data analysis for something like 10 years now (gulp!) and teaching data analysis in person for 6+ years. One of the things we do in my data analysis class at Hopkins is to perform a complete data analysis (from raw data to written report) every couple of weeks. Then I grade each assignment for everything from data cleaning to the written report and reproducibility. I’ve noticed over the course of teaching this class (and classes online) that there are many common elements of data analytic style that I don’t often see in textbooks, or when I do, I see them spread across multiple books.

I’ve posted on some of these issues in some open source guides I’ve posted to Github like:

• 10 things statistics taught us about big data analysis
• The Leek Group Guide to R packages
• How to share data with a statistician

But I decided that it might be useful to have a more complete guide to the “art” part of data analysis. One goal is to summarize in a succinct way the most common difficulties encountered by practicing data analysts. It may be a useful guide for peer reviewers who could refer to section numbers when evaluating manuscripts, for instructors who have to grade data analyses, as a supplementary text for a data analysis class, or just as a useful reference. It is modeled loosely in format and aim on the Elements of Style by William Strunk. Just as with the EoS, both the checklist and my book cover a small fraction of the field of data analysis, but my experience is that once these elements are mastered, data analysts benefit most from hands on experience in their own discipline of application, and that many principles may be non-transferable beyond the basics. But just as with writing, new analysts would do better to follow the rules until they know them well enough to violate them.

• Buy EDAS on Leanpub
• Buy EDAS on Amazon

The book includes a basic checklist that may be useful as a guide for beginning data analysts or as a rubric for evaluating data analyses. I’m reproducing it here so you can comment/hate/enjoy on it.

The data analysis checklist

This checklist provides a condensed look at the information in this book. It can be used as a guide during the process of a data analysis, as a rubric for grading data analysis projects, or as a way to evaluate the quality of a reported data analysis.

I Answering the question

1. Did you specify the type of data analytic question (e.g. exploration, assocation causality) before touching the data?

2. Did you define the metric for success before beginning?

3. Did you understand the context for the question and the scientific or business application?

4. Did you record the experimental design?

5. Did you consider whether the question could be answered with the available data?

II Checking the data

1. Did you plot univariate and multivariate summaries of the data?

2. Did you check for outliers?

3. Did you identify the missing data code?

III Tidying the data

1. Is each variable one column?

2. Is each observation one row?

3. Do different data types appear in each table?

4. Did you record the recipe for moving from raw to tidy data?

5. Did you create a code book?

6. Did you record all parameters, units, and functions applied to the data?

IV Exploratory analysis

1. Did you identify missing values?

2. Did you make univariate plots (histograms, density plots, boxplots)?

3. Did you consider correlations between variables (scatterplots)?

4. Did you check the units of all data points to make sure they are in the right range?

5. Did you try to identify any errors or miscoding of variables?

6. Did you consider plotting on a log scale?

7. Would a scatterplot be more informative?

V Inference

1. Did you identify what large population you are trying to describe?

2. Did you clearly identify the quantities of interest in your model?

3. Did you consider potential confounders?

4. Did you identify and model potential sources of correlation such as measurements over time or space?

5. Did you calculate a measure of uncertainty for each estimate on the scientific scale?

VI Prediction

1. Did you identify in advance your error measure?

2. Did you immediately split your data into training and validation?

3. Did you use cross validation, resampling, or bootstrapping only on the training data?

4. Did you create features using only the training data?

5. Did you estimate parameters only on the training data?

6. Did you fix all features, parameters, and models before applying to the validation data?

7. Did you apply only one final model to the validation data and report the error rate?

VII Causality

1. Did you identify whether your study was randomized?

2. Did you identify potential reasons that causality may not be appropriate such as confounders, missing data, non-ignorable dropout, or unblinded experiments?

3. If not, did you avoid using language that would imply cause and effect?

VIII Written analyses

1. Did you describe the question of interest?

2. Did you describe the data set, experimental design, and question you are answering?

3. Did you specify the type of data analytic question you are answering?

4. Did you specify in clear notation the exact model you are fitting?

5. Did you explain on the scale of interest what each estimate and measure of uncertainty means?

6. Did you report a measure of uncertainty for each estimate on the scientific scale?

IX Figures

1. Does each figure communicate an important piece of information or address a question of interest?

2. Do all your figures include plain language axis labels?

3. Is the font size large enough to read?

4. Does every figure have a detailed caption that explains all axes, legends, and trends in the figure?

X Presentations

1. Did you lead with a brief, understandable to everyone statement of your problem?

2. Did you explain the data, measurement technology, and experimental design before you explained your model?

3. Did you explain the features you will use to model data before you explain the model?

4. Did you make sure all legends and axes were legible from the back of the room?

XI Reproducibility

1. Did you avoid doing calculations manually?

2. Did you create a script that reproduces all your analyses?

3. Did you save the raw and processed versions of your data?

4. Did you record all versions of the software you used to process the data?

5. Did you try to have someone else run your analysis code to confirm they got the same answers?

XI R packages

1. Did you make your package name “Googleable”

2. Did you write unit tests for your functions?

3. Did you write help files for all functions?

4. Did you write a vignette?

5. Did you try to reduce dependencies to actively maintained packages?

6. Have you eliminated all errors and warnings from R CMD CHECK?

Advanced Statistics for the Life Sciences MOOC Launches Today

In In we will teach statistical techniques that are commonly used in the analysis of high-throughput data and their corresponding R implementations. In Week 1 we will explain inference in the context of high-throughput data and introduce the concept of error controlling procedures. We will describe the strengths and weakness of the Bonferroni correction, FDR and q-values. We will show how to implement these in cases in which  thousands of tests are conducted, as is typically done with genomics data. In Week 2 we will introduce the concept of mathematical distance and how it is used in exploratory data analysis, clustering, and machine learning. We will describe how techniques such as principal component analysis (PCA) and the singular value decomposition (SVD) can be used for dimension reduction in high dimensional data. During week 3 we will describe confounding, latent variables and factor analysis in the context of high dimensional data and how this relates to batch effects. We will show how to implement methods such as SVA to perform inference on data affected by batch effects. Finally, during week 4 we will show how statistical modeling, and empirical Bayes modeling in particular, are powerful techniques that greatly improve precision in high-throughput data. We will be using R code to explain concepts throughout the course. We will also be using exploratory data analysis and data visualization to motivate the techniques we teach during each week.

Navigating Big Data Careers with a Statistics PhD