ThinkStats … in R (including Example 1.2)

March 4, 2012
By

(This article was first published on rud.is » R, and kindly contributed to R-bloggers)

ThinkStats (by Allen B. Downey) is a good book to get you familiar with statistics (and even Python, if you’ve done some scripting in other languages).

I thought it would be interesting to present some of the examples & exercises in the book in R. Why? Well, once you’ve gone through the material in a particular chapter the “hard way”, seeing how you’d do the same thing in a language specifically designed for statistical computing should show when it’s best to use such a domain specific language and when you might want to consider a hybrid approach. I am also hoping it helps make R a bit more accessible to folks.

You’ll still need the book and should work through the Python examples to get the most out of these posts.

I’ll try to get at least one example/exercise section up a week.

Please submit all errors, omissions or optimizations in the comments section.

The star of the show is going to be the “data frame” in most of the examples (and is in this one). Unlike the Python code in the book, most of the hard work here is figuring out how to use the data frame file reader to parse the ugly fields in the CDC data file. By using some tricks, we can approximate the “field start:length” style of the Python code but still keep the automatic reading/parsing of the R code (including implicit handling of “NA” values).

The power & simplicity of using R’s inherent ability to apply a calculation across a whole column (pregnancies$agepreg <- pregnancies$agepreg / 100) should jump out. Unfortunately, not all elements of the examples in R will be as nice or straightforward.

You'll also notice that I cheat and use str() for displaying summary data.

Enough explanation! Here's the code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
# ThinkStats in R by @hrbrmstr
# Example 1.2
# File format info: http://www.cdc.gov/nchs/nsfg/nsfg_cycle6.htm
 
# setup a data frame that has the field start/end info
 
pFields <- data.frame(name  = c('caseid', 'nbrnaliv', 'babysex', 'birthwgt_lb','birthwgt_oz','prglength', 'outcome', 'birthord',  'agepreg',  'finalwgt'), 
                      begin = c(1, 22, 56, 57, 59, 275, 277, 278, 284, 423), 
                      end   = c(12, 22, 56, 58, 60, 276, 277, 279, 287, 440) 
) 
 
# calculate widths so we can pass them to read.fwf()
 
pFields$width <- pFields$end - pFields$begin + 1 
 
# we aren't reading every field (for the book exercises)
 
pFields$skip <-  (-c(pFields$begin[-1]-pFields$end[-nrow(pFields)]-1,0)) 
 
widths <- c(t(pFields[,4:5])) 
widths <- widths[widths!=0] 
 
# read in the file
 
pregnancies <- read.fwf("2002FemPreg.dat", widths) 
 
# assign column names
 
names(pregnancies) <- pFields$name 
 
# divide mother's age by 100
 
pregnancies$agepreg <- pregnancies$agepreg / 100
 
# convert weight at birth from lbs/oz to total ounces
 
pregnancies$totalwgt_oz = pregnancies$birthwgt_lb * 16 + pregnancies$birthwgt_oz
 
rFields <- data.frame(name  = c('caseid'), 
                      begin = c(1), 
                      end   = c(12) 
) 
 
rFields$width <- rFields$end - rFields$begin + 1 
rFields$skip <-  (-c(rFields$begin[-1]-rFields$end[-nrow(rFields)]-1,0)) 
 
widths <- c(t(rFields[,4:5])) 
widths <- widths[widths!=0] 
 
respondents <- read.fwf("2002FemResp.dat", widths) 
names(respondents) <- rFields$name
 
str(respondents)
str(pregnancies)

To leave a comment for the author, please follow the link and comment on his blog: rud.is » R.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...



If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Tags: , , , , , , , , , ,

Comments are closed.