Wednesday, April 30, 2014

Three Ways to Format R Code for Blogger

Unformatted Code

If you are like me originally then you might not think it is worth it to spend the extra energy to format your code.  After all people can just copy what you have and paste it into their preferred editor which will do its own formatting, no sweat.

Well, this might be true for some code in some languages, but in R it really is not gong to fly.  The problem is the ubiquitous <- which is very easily confused with html tags which take the form <tag>. Thus, you are going to be looking for a way to make your code not be buggy when you post it.  One way might be to find and replace all < with &lt; and > with &gt; using a text editor but that is already a fair amount of trouble especially if you actually want to have some HTML tags in your post.  Besides, if you are going to all of that trouble, you might as well as well actually format your posts with syntax highlighting since it will be that much more enjoyable for readers.

Here are three simple ways:

1. Export Code through Notepad++
Exporting through Notepad++ is easy and embeds a css style sheet within the HTML document which makes it easy to edit styles after export.  There is also a fair amount of editing possible before export as well.
  1. Open myscript.R file
  2. Configure appearance of formatting with: Settings>Style>R to your taste
  3. Export HTML: Plugins>NPPExport>Export to HTML
  4. Open saved myscript.html as a text file.
  5. Copy and paste HTML into new post HTML in blogger. 

2. Format with Pretty-R
Pretty-R is useful because it is very simple to use. It also provides links from highlighted function to associated R documentation.  However, the style is fixed and difficult to modify before or after export.
  1. Go to
  2. Paste your R script
  3. Copy HTML code in the HTML box and paste into new post HTML in blogger.

3. Format with knittr
knittr is a pretty amazing package which turns R markdown files into HTML files.  There is a lot of customization possible with it and it is under active development.  I have only recently started looking into it though it seems to me to be the way to go forward. In addition, RStudio supports knittr markdown formatting.  However, it does take a few extra steps to set up but ultimately gives you more control over your final output.
  1. Open RStudio
  2. Select File>New>R Markdown
  3. Paste your code between ```{r} and ``` and structure appropriately.
  4. Specify formatting of your posts as you see fit. There are a lot of great options.  Images will be embedded directly into the HTML document.  This is a pretty amazing feature since it means your file will not require any external dependencies.
  5. Click Knit HTML. This should create an HTML document identical of the same name as the Rmd file.
  6. Edit the generated HTML file.
  7. You may need to remove the ccs classes h1 through h6 as they can interfere with your blogger theme's css classes of the same name.  You may also need to remove the markdown titled identified as <h1> (around line 150-200) if you are going to input your post title in the default blogger title position (which is recommended).
  8. Copy you edited HTML code into blogger

(Note: I know step 7 is a bit hacky so I have asked the question on how to change the default css styles with knitr so that they are not interfering with blogger on stackoverflow.  Please take a look at the question as I hope someone will be able to answer it soon and greatly simply this process.)
Flattr this

Tuesday, April 29, 2014

Dave Giles on "MCMC for Econometrics Students"

In an excellent four part series of posts in March, Dave Giles introduces Markov Chain Monte Carlo (MCMC) and Gibbs samplers.  In these posts he gives a cogent explanation for the reasoning and mechanics involved in this branch of econometrics/statistics as well as clear simulated examples in R.

If you have not checked it out yet, now is definitely a good time.

Find Dave's posts:
Post 1: Introduction
Post 2: Showing that MCMC "Actually Works"
Post 3: Shows an additional example as well as how to extract marginal posterior distributions
Post 4: Shows how simple it is to use R to implement MCMC

In order to experiment with the topic of MCMC I have made some modifications to Dave's code in R.  He makes no assertions that his code is in efficient form.

Gibs defined below is the same as his code except that Gibs is now defined as a function.  Gibs2 has modified the code as best I could do in order so that I am working with vectors as much as possible rather than item by item manipulation.  I used Noam Ross's excellent post to inform out my understanding of improving processing speeds with R.

By vectoring the random draws Gibs2 processes 2 to 3 times faster than Gibs.  Full code can be found on github:

# user system elapsed
# 0.97 0.06   1.17

# user system elapsed
# 9.31 0.02   9.64

# user system elapsed
# 0.35 0.00   0.36
# user system elapsed
# 3.66 0.01   3.80 

As an exercise I also coded the same gibs function in Julia.  This can also be found on github as well.

@elapsed yy = gibs()
# 0.063151467
@elapsed gibs(10^6)
# 0.479542057

@elapsed yy = gibs2()
# 0.010729382
@elapsed yy = gibs2(10^6)
# 0.065821774
The first thing to notice is that when coding the initial form of gibs, the julia version is considerably faster (10-15x).   Gains from improving the form are also larger in julia with gibs2 running much faster than the R version (30-55x faster).
Flattr this

Monday, April 28, 2014

Tall Stata, low math, extra pixels

With a clever titled blog and in well written in-depth initial post, Alex Gamma at the University of Zurich, Switzerland enters the Stata blogging scene. Check out his post on "The nicest place to live in Switzerland" at
Tall Stata
 Find a list of syndicated Stata sites at

Thursday, April 24, 2014

Shout out to "R Handles Big Data"

Searching for Significance

I just wanted to give a shout out to check out this post by Bob Muenchen on his excellent blog for his exceptional post a little over a year ago entitle "R Handles Big Data".

Wednesday, April 23, 2014

A Weekend With Julia: An R User's Reflections
The Famous Julia

First off, I am not going to talk much about Julia's speed. Everybody has seen the tables and graphs showing how in this benchmark or another, Julia is tens times or a hundred times faster than R.  Most blog posts talking about Julia test the generality of these results (Bogumił Kamiński 2013, Randy Zwitch 2013, and  Wes McKinney 2012).

Enough said about machine speed!  Let's talk more about intuitive appeal, compactness of notation, and aesthetics.

Julia has some very thoughtful design features that make it an extremely enjoyable language to program in despite being in a nascent state.

1. Julia has some great ways of handling string expressions.  In Julia all strings are subsettable.  Thus:

julia> a = "Hello world"; a[3:7]
"llo w"

In R:
R>a <- "Hello world"; substr(a, 3, 7)

Also, "Julia allows interpolation into string literals using $, as in Perl:"(doc)

julia>user = "Bob"
julia>"Hello $user. How are you?"
"Hello Bob. How are you?"

2. Julia implements comprehension syntax for defining arrays which is an incredibly powerful method. Formally it looks something like this: A = [ F(x,y,...) for x=rx, y=ry, ... ]

Imagine we would like to define an area equivalent to the number line:

julia> A = [ x*y for x=1:12, y=1:12]; A
12x12 Array{Int64,2}:
  1   2   3   4   5   6   7   8    9   10   11   12
  2   4   6   8  10  12  14  16   18   20   22   24
  3   6   9  12  15  18  21  24   27   30   33   36
  4   8  12  16  20  24  28  32   36   40   44   48
  5  10  15  20  25  30  35  40   45   50   55   60
  6  12  18  24  30  36  42  48   54   60   66   72
  7  14  21  28  35  42  49  56   63   70   77   84
  8  16  24  32  40  48  56  64   72   80   88   96
  9  18  27  36  45  54  63  72   81   90   99  108
 10  20  30  40  50  60  70  80   90  100  110  120
 11  22  33  44  55  66  77  88   99  110  121  132
 12  24  36  48  60  72  84  96  108  120  132  144

Let's see how we might do this in R.
R><-matrix(nrow=12,ncol=12); A <-col(A)*row(A); A
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
 [1,]    1    2    3    4    5    6    7    8    9    10    11    12
 [2,]    2    4    6    8   10   12   14   16   18    20    22    24
 [3,]    3    6    9   12   15   18   21   24   27    30    33    36
 [4,]    4    8   12   16   20   24   28   32   36    40    44    48
 [5,]    5   10   15   20   25   30   35   40   45    50    55    60
 [6,]    6   12   18   24   30   36   42   48   54    60    66    72
 [7,]    7   14   21   28   35   42   49   56   63    70    77    84
 [8,]    8   16   24   32   40   48   56   64   72    80    88    96
 [9,]    9   18   27   36   45   54   63   72   81    90    99   108
[10,]   10   20   30   40   50   60   70   80   90   100   110   120
[11,]   11   22   33   44   55   66   77   88   99   110   121   132
[12,]   12   24   36   48   60   72   84   96  108   120   132   144

3. Functions can be written in a mathematically intuitive form.
julia> f(x,y)=x^3-y+x*y; f(3,2)

4. Numerical constants leading functions are automatically interpreted
julia> x=13; 3x

julia> (12+4)x

R> x<-13; 3*x

R> (12+4)*x

5. Julia does not worry about deep nesting of functions. Imagine a very silly function that adds up all of the integers between 1 and n.
julia> f(n)=(n>1 ? f(n-1)+n: n); f(100)

R> f <-function(n) ifelse(n>1, f(n-1)+n, n); f(100)
[1] 5050

Both R and Julia seem to work fine.  But what happens when we go a little deeper?
R> f(10^3)
Fails while

julia> f(10^5)

Does not even cause a hiccup! Now I have never been in the situation of needing this many recursions but this result reflects the general power of the language.

6. Julia has no problem with many Unicode characters. Thus if you want θ(μ,σ,φ)=μ^σ/φ
julia> θ(μ,σ,φ)=μ^σ/φ; θ(1,2,3)

Personally I find this notation extremely appealing as it succinctly communicates the equation for which the researcher is actually dealing with as opposed to what typically happens when I code R.

R> theta <function(mu, sigma, phi) = mu^sigma/phi; theta(1,2,3)

Besides being more compact, Julia's notation requires less mental juggling in order to translate from mathematical concepts to direct implementation.

7. Tuple fast assignment.

I am not sure if I am referring to this properly. Perhaps one of the omniscient badasses working on Julia can correct me if I am not describing this feature correctly.

I will just show you how it works:
julia> a, b, c, d = 1, 22, 33, 444
julia> b,d
(22, 444)
This might seem very trivial but check out how easy it is to pass parameters into a function. Item response theory probability of a 4PL model defined with ability parameter θ and item parameters a,b,c,d.

function ppi(θ=0, ipar=(1, 0, 0, 1), D=1.7)
  a, b, c, d = ipar
  # Define exponential factor appearing in
  # numerator and denominator.
  ε = e^(D * a *  - b))
  # Define final probability being returned.
  c + ((d-c)*ε)/(1+ε)

julia> ppi()

8. The Readability of programming in Julia.

What I mean by this is that when programming in Julia, so long as you follow basic rules regarding programming efficiency your program is likely to run at a comparable speed in Julia as it would be if it were programmed in C. For me this is a big deal. It means that I do not have to think about figuring out how to program the same task in two different languages if I am programming a new library or feature.

Likewise, when reading functions and libraries written in Julia, I hopefully will have to worry less about working with other languages. This should make it possible for the source code for functions to be more accessible for the purposes of learning and improving my understanding of programming.

Some Critiques
It is at this time that I want to mention a few critiques.  As I see it, these critiques are entirely based on the novelty of the Julia environment which is completely addressed by the creators listing the current version of the language in 0.2 and testing 0.3.  As far as how the language actually performs, I have no complaints.  It does not have as many libraries as R but they are in development.

1. The documentation is surprisingly very readable and surprisingly enjoyable.  I found that reading through documentation not only enhanced my understanding of Julia but significantly enhanced my understanding of programming languages in general.  That said, the documentation seems to be written from programming polyglots to other programming polyglots.  I think it would be very difficult for someone without programming experience to read the documentation and be able to do much with Julia.

2. The documentation is sparse on examples and the examples tend to be at a high level of complexity.  I personally think that the best way to learn to program is through examples (though this might explain my spotting understanding of abstract concepts like scope).  The documentation of Julia has a fair number of examples but I would think that doubling the number of examples would be very helpful.  In addition, there are quite a number of functions written for Julia which do not have any documentation at all.  Obviously it would be helpful to have documentation for such functions.

3. Finally, as a windows users (I know, poor form) it is not very easy to use Julia.  I have been using Forio's Julia Studio (since I do not like working with Window's command line) but don't find it that great to work with plus it is proprietary (though free).  It would be nice if there was a very basic windows IDE like R's for which I could write code in Nopepad++ and send it to Julia.  But once again this is a minor issue.  I am resolved to get a Linux build running on my machine once I make my way back to a country with decent internet.