R Basics for Lecture 1
Econ B2000, Statistics and Introduction to Econometrics
Kevin R Foster, the City College of New York, CUNY
For this class we’ll be using the statistical analysis program, R. The computer labs here on campus have the necessary software. If you have your own computer, download and install the program from R Project along with R-Studio. Depending on your particular device these details will be different so take it patiently and figure it out - there are sufficient online sources for help, to address any problem that arises.
The best way to learn it is to just do it. How many times have you got a new game, skipped reading the instructions, and learned by crashing it a few times? This isn’t quite so simple but the basic gist remains - try it. Below I give you some pointers about how to just do it, but they’re meant to be read then immediately done in real life, I’ll give you some commands to just copy-and-paste into the program.
The program, R, is the machine underneath that does the work, while R-Studio is a skin on top that makes it easier to use. Install both then run R-Studio, and you’ll get something that looks like this: 
The screen has 4 parts: the Console at lower left (where I drew the green arrow) is most important since that’s where you type in the commands to R.
You can start by just copying and pasting commands from this help into the “Console” and seeing the output from the program.
The guide, An Introduction to R, suggests these commands to give a basic flavor of what’s done and hint at the power. Don’t worry too much about each step for now, this is just a glimpse to show that with a few commands you can do some relatively sophisticated estimation – i.e. for a small cost you can get a large benefit. Copy them and paste into the “Console”.
x <- 1:50
w <- 1 + sqrt(x)/2
example1 <- data.frame(x=x, y= x + rnorm(x)*w)
attach(example1)This creates x and y variables (where the rnorm command creates random numbers from a normal distribution), puts them into a data frame, then attaches that data frame so that R can use it in later calculations. (What’s a data frame? A way to collect together related data items.) Next some stats - create a linear model (that’s the “lm”) then a “lowess” nonparametric local regression, and plot the two estimations to compare. (Just copy and paste these, don’t worry about understanding them for now!)
fm <- lm(y ~ x)
summary(fm)##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.9479 -2.0413 0.5088 2.4289 5.3096
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.32772 0.95071 0.345 0.732
## x 0.98486 0.03245 30.353 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.311 on 48 degrees of freedom
## Multiple R-squared: 0.9505, Adjusted R-squared: 0.9494
## F-statistic: 921.3 on 1 and 48 DF, p-value: < 2.2e-16lrf <- lowess(x, y)
plot(x, y)
lines(x, lrf$y)
abline(0, 1, lty=3)
abline(coef(fm))
detach()Your own graph should look similar although with the random numbers, not exactly.
The final “detach” command just cleans up, it is the opposite of “attach”.
Later you will find that while it is possible to go back and fix up code that you previously ran in R (the up arrow brings most recent commands), it is a bit of a pain. It is easier to write out the code in the upper-left panel, then R-Studio can obediently run all or part of it in R for you. Then if you make a mistake or just want to do it again, it’s a bit easier. You can save your code for next time, as well.
One big thing to learn is that while pushing buttons is easier at first, eventually you want to be able to write code. I won’t insist that you do that from week #1, but keep it in mind as a goal to work towards. In R, you can usually see the code generated by the button-pushes so you can learn.
For all of these commands, you can use the help box on R Studio. But as I said, don’t worry much about those commands for now, I’m not expecting you to become an R-ninja overnight.
With Some Data
Next we will go through some basic stats with R using the Census Bureau’s PUMS (Public Use Microdata Sample, from the American Community Survey, accessed from IPUMS). I have restricted the sample to contain only people living in the state of New York.
You have to learn a bit about your computer’s filing system. When you download stuff it probably goes into a Download folder. When you installed R and R-Studio, those went into their own folders (probably within some Application folder) and then the program might have created a R folder for your work. When you do homework or other projects, those might (well, really should) get their own folder. But joining those up is a bit of housekeeping.
Go and download the PUMS data from the class page, that will put that zip file into your Downloads folder. Within that zip file is one particular file, pums_NY.RData - move that into some appropriately-titled folder for this project. (Personally I might use the default “R Projects” folder then create “PUMS Data” as a folder within that, but you can choose. Young people with supple minds can get by even if they have zero organization like a hippie, but later you might appreciate being tidier.)
R is looking for files in a particular directory. Type the command, “getwd()” to see where it’s currently looking. Then use the command, “setwd,” to tell it where it ought to be looking (the R Projects/PUMS Data folder that you created just now, or wherever you put that data file). Alt you can click “Session” then “Set Working Directory” then “Choose Directory” and click the folder, and that will insert the “setwd” command onto the Console line. Tell R to look for the data in the folder where you put the data.
If you have a Dropbox account then that can be worked in the same way, if you set that up as a folder on your computer.
Then run these commands (output from those commands is below),
rm(list = ls(all = TRUE)) # clear workspace
load("pums_NY.RData")
str(dat_pums_NY)## 'data.frame': 196314 obs. of 57 variables:
## $ Age : num 43 45 33 57 52 26 83 87 21 45 ...
## $ female : num 1 0 0 0 1 0 1 0 0 0 ...
## $ PERNUM : num 1 2 1 1 2 3 1 2 1 1 ...
## $ educ_nohs : num 0 0 0 0 1 0 0 0 0 0 ...
## $ educ_hs : num 0 0 1 1 0 0 1 1 0 0 ...
## $ educ_smcoll : num 0 1 0 0 0 0 0 0 1 0 ...
## $ educ_as : num 0 0 0 0 0 1 0 0 0 1 ...
## $ educ_bach : num 0 0 0 0 0 0 0 0 0 0 ...
## $ educ_adv : num 1 0 0 0 0 0 0 0 0 0 ...
## $ ANCESTR1D : num 2610 511 880 7060 7060 7060 1440 1420 9020 321 ...
## $ ANCESTR2D : num 9990 9990 9990 9990 9990 9990 9990 9990 9990 9990 ...
## $ immig : num 0 0 0 1 1 1 0 0 0 0 ...
## $ Hispanic : num 1 0 0 0 0 0 0 0 0 0 ...
## $ Hisp_Mex : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Hisp_PR : num 1 0 0 0 0 0 0 0 0 0 ...
## $ Hisp_Cuban : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Hisp_DomR : num 0 0 0 0 0 0 0 0 0 0 ...
## $ white : num 1 1 1 0 0 0 1 1 0 1 ...
## $ AfAm : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Amindian : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Asian : num 0 0 0 1 1 1 0 0 0 0 ...
## $ race_oth : num 0 0 0 0 0 0 0 0 1 0 ...
## $ Married : num 1 1 0 1 1 0 1 1 0 1 ...
## $ divwidsep : num 0 0 0 0 0 0 0 0 0 0 ...
## $ unmarried : num 0 0 1 0 0 1 0 0 1 0 ...
## $ veteran : num 0 0 0 0 0 0 0 0 0 0 ...
## $ has_AnyHealthIns : num 1 1 1 0 0 0 1 1 1 1 ...
## $ has_PvtHealthIns : num 1 1 1 0 0 0 1 1 1 1 ...
## $ Commute_car : num 1 1 1 0 0 0 0 0 1 1 ...
## $ Commute_bus : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Commute_subway : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Commute_rail : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Commute_other : num 0 0 0 0 1 0 0 0 0 0 ...
## $ below_povertyline: num 0 0 0 1 1 1 0 0 0 0 ...
## $ below_150poverty : num 0 0 0 1 1 1 0 0 0 0 ...
## $ below_200poverty : num 0 0 0 1 1 1 1 1 0 0 ...
## $ foodstamps : num 1 1 1 1 1 1 1 1 1 1 ...
## $ work_fullyr : num 1 1 1 0 0 0 0 0 1 1 ...
## $ income_total : num 110000 39000 72000 0 7000 0 5800 19200 5200 61000 ...
## $ income_wagesal : num 110000 39000 72000 0 7000 0 0 0 5200 61000 ...
## $ HH_income : num 0 0 72000 7000 7000 7000 25000 25000 0 81000 ...
## $ owner_cost : num 2850 2850 0 0 0 ...
## $ rent_cost : num 0 0 430 900 900 900 0 0 0 0 ...
## $ occ_dum : num 1820 1550 4210 0 4520 0 0 0 4760 6440 ...
## $ ind_dum : num 7860 3390 770 0 8980 0 0 0 5170 770 ...
## $ in_NYC : num 0 0 0 1 1 1 0 0 0 0 ...
## $ PUMA : num 3106 3106 100 4103 4103 ...
## $ in_Bronx : num 0 0 0 0 0 0 0 0 0 0 ...
## $ in_Manhattan : num 0 0 0 0 0 0 0 0 0 0 ...
## $ in_StatenI : num 0 0 0 0 0 0 0 0 0 0 ...
## $ in_Brooklyn : num 0 0 0 0 0 0 0 0 0 0 ...
## $ in_Queens : num 0 0 0 1 1 1 0 0 0 0 ...
## $ in_Westchester : num 1 1 0 0 0 0 0 0 0 0 ...
## $ in_Nassau : num 0 0 0 0 0 0 0 0 0 0 ...
## $ ROOMS : num 8 8 2 3 3 3 7 7 0 5 ...
## $ BUILTYR2 : num 10 10 9 5 5 5 4 4 0 6 ...
## $ UNITSSTR : num 4 4 3 10 10 10 3 3 0 3 ...attach(dat_pums_NY)The command, str, is a good way to check out the data. In the next section I’ll explain more about the data and what the lines mean but for now Age is the person’s age in years and female is a logical 0/1 variable. So if you look at the output from str, see that the first person is a 43-year-old female, next is a 45-year-old male, etc.
str(Age)## num [1:196314] 43 45 33 57 52 26 83 87 21 45 ...str(female)## num [1:196314] 1 0 0 0 1 0 1 0 0 0 ...You can also use the command, summary, to find out about data.
summary(dat_pums_NY)## Age female PERNUM educ_nohs
## Min. : 0.00 Min. :0.0000 Min. : 1.000 Min. :0.0000
## 1st Qu.:20.00 1st Qu.:0.0000 1st Qu.: 1.000 1st Qu.:0.0000
## Median :41.00 Median :1.0000 Median : 2.000 Median :0.0000
## Mean :40.58 Mean :0.5191 Mean : 2.122 Mean :0.3166
## 3rd Qu.:59.00 3rd Qu.:1.0000 3rd Qu.: 3.000 3rd Qu.:1.0000
## Max. :94.00 Max. :1.0000 Max. :16.000 Max. :1.0000
## educ_hs educ_smcoll educ_as educ_bach
## Min. :0.0000 Min. :0.0000 Min. :0.00000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.0000
## Median :0.0000 Median :0.0000 Median :0.00000 Median :0.0000
## Mean :0.2196 Mean :0.1576 Mean :0.06641 Mean :0.1374
## 3rd Qu.:0.0000 3rd Qu.:0.0000 3rd Qu.:0.00000 3rd Qu.:0.0000
## Max. :1.0000 Max. :1.0000 Max. :1.00000 Max. :1.0000
## educ_adv ANCESTR1D ANCESTR2D immig
## Min. :0.0000 Min. : 10 Min. : 10 Min. :0.00000
## 1st Qu.:0.0000 1st Qu.: 511 1st Qu.:6801 1st Qu.:0.00000
## Median :0.0000 Median :2251 Median :9990 Median :0.00000
## Mean :0.1025 Mean :3952 Mean :7668 Mean :0.08103
## 3rd Qu.:0.0000 3rd Qu.:9020 3rd Qu.:9990 3rd Qu.:0.00000
## Max. :1.0000 Max. :9990 Max. :9990 Max. :1.00000
## Hispanic Hisp_Mex Hisp_PR Hisp_Cuban
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. :0.000000
## 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.000000
## Median :0.0000 Median :0.00000 Median :0.0000 Median :0.000000
## Mean :0.1353 Mean :0.01672 Mean :0.0454 Mean :0.003107
## 3rd Qu.:0.0000 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:0.000000
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :1.000000
## Hisp_DomR white AfAm Amindian
## Min. :0.00000 Min. :0.0000 Min. :0.0000 Min. :0.000000
## 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.000000
## Median :0.00000 Median :1.0000 Median :0.0000 Median :0.000000
## Mean :0.02642 Mean :0.7093 Mean :0.1347 Mean :0.003994
## 3rd Qu.:0.00000 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:0.000000
## Max. :1.00000 Max. :1.0000 Max. :1.0000 Max. :1.000000
## Asian race_oth Married divwidsep
## Min. :0.00000 Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :0.00000 Median :0.0000 Median :0.0000 Median :0.0000
## Mean :0.07088 Mean :0.1209 Mean :0.3937 Mean :0.1515
## 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:0.0000
## Max. :1.00000 Max. :1.0000 Max. :1.0000 Max. :1.0000
## unmarried veteran has_AnyHealthIns has_PvtHealthIns
## Min. :0.0000 Min. :0.00000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:1.0000 1st Qu.:0.0000
## Median :0.0000 Median :0.00000 Median :1.0000 Median :1.0000
## Mean :0.4548 Mean :0.05515 Mean :0.9082 Mean :0.6706
## 3rd Qu.:1.0000 3rd Qu.:0.00000 3rd Qu.:1.0000 3rd Qu.:1.0000
## Max. :1.0000 Max. :1.00000 Max. :1.0000 Max. :1.0000
## Commute_car Commute_bus Commute_subway Commute_rail
## Min. :0.0000 Min. :0.00000 Min. :0.00000 Min. :0.00000
## 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.00000 1st Qu.:0.00000
## Median :0.0000 Median :0.00000 Median :0.00000 Median :0.00000
## Mean :0.2892 Mean :0.02343 Mean :0.06291 Mean :0.01165
## 3rd Qu.:1.0000 3rd Qu.:0.00000 3rd Qu.:0.00000 3rd Qu.:0.00000
## Max. :1.0000 Max. :1.00000 Max. :1.00000 Max. :1.00000
## Commute_other below_povertyline below_150poverty below_200poverty
## Min. :0.00000 Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :0.00000 Median :0.0000 Median :0.0000 Median :0.0000
## Mean :0.05029 Mean :0.1434 Mean :0.2238 Mean :0.3012
## 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:0.0000 3rd Qu.:1.0000
## Max. :1.00000 Max. :1.0000 Max. :1.0000 Max. :1.0000
## foodstamps work_fullyr income_total income_wagesal
## Min. :0.0000 Min. :0.000 Min. : -11300 Min. : 0
## 1st Qu.:1.0000 1st Qu.:0.000 1st Qu.: 0 1st Qu.: 0
## Median :1.0000 Median :0.000 Median : 14000 Median : 0
## Mean :0.8351 Mean :0.376 Mean : 31033 Mean : 22646
## 3rd Qu.:1.0000 3rd Qu.:1.000 3rd Qu.: 40000 3rd Qu.: 30000
## Max. :1.0000 Max. :1.000 Max. :1125000 Max. :504000
## HH_income owner_cost rent_cost occ_dum
## Min. :-12800 Min. : 0 Min. : 0.0 Min. : 0
## 1st Qu.: 0 1st Qu.: 0 1st Qu.: 0.0 1st Qu.: 0
## Median : 21400 Median : 702 Median : 0.0 Median :1030
## Mean : 30061 Mean : 1148 Mean : 336.2 Mean :2474
## 3rd Qu.: 55000 3rd Qu.: 1825 3rd Qu.: 550.0 3rd Qu.:4700
## Max. : 99990 Max. :12248 Max. :3700.0 Max. :9920
## ind_dum in_NYC PUMA in_Bronx
## Min. : 0 Min. :0.0000 Min. : 100 Min. :0.00000
## 1st Qu.: 0 1st Qu.:0.0000 1st Qu.:1500 1st Qu.:0.00000
## Median :3980 Median :0.0000 Median :3107 Median :0.00000
## Mean :3893 Mean :0.3571 Mean :2701 Mean :0.05582
## 3rd Qu.:7860 3rd Qu.:1.0000 3rd Qu.:3901 3rd Qu.:0.00000
## Max. :9920 Max. :1.0000 Max. :4114 Max. :1.00000
## in_Manhattan in_StatenI in_Brooklyn in_Queens
## Min. :0.00000 Min. :0.00000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.00000 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :0.00000 Median :0.00000 Median :0.0000 Median :0.0000
## Mean :0.04695 Mean :0.02036 Mean :0.1254 Mean :0.1086
## 3rd Qu.:0.00000 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:0.0000
## Max. :1.00000 Max. :1.00000 Max. :1.0000 Max. :1.0000
## in_Westchester in_Nassau ROOMS BUILTYR2
## Min. :0.00000 Min. :0.0000 Min. : 0.000 Min. : 0.000
## 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.: 4.000 1st Qu.: 1.000
## Median :0.00000 Median :0.0000 Median : 6.000 Median : 3.000
## Mean :0.04385 Mean :0.0691 Mean : 5.825 Mean : 3.442
## 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.: 7.000 3rd Qu.: 5.000
## Max. :1.00000 Max. :1.0000 Max. :19.000 Max. :17.000
## UNITSSTR
## Min. : 0.00
## 1st Qu.: 3.00
## Median : 3.00
## Mean : 4.34
## 3rd Qu.: 6.00
## Max. :10.00print(NN_obs <- length(Age))## [1] 196314So this shows that there are 196,314 people in this dataset.
Simple Stats
We compare the average age of the men and the women in the data,
summary(Age[female == 1])## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 21.00 43.00 41.88 60.00 94.00summary(Age[!female])## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 19.00 39.00 39.19 57.00 94.00The female dummy variable is a logical term, with a zero or one for false or true. The comparison is between those who have the variable female=1 (i.e. women) and those not female=1 (so logical not, denoted with the “!” symbol, i.e. men). I know, you can – and people do! – worry that this binary classification for gender misses some people; government statistics are just not there yet.
So women are, on average, a bit older, with an average age of 41.9 compared with 39.2 for men. You might wonder (if you were to begin to think like a statistician) whether that is a big difference - hold onto that thought! Alternately you can use the commands to calculate the average, mean(), and the standard deviation, sd(), to get those statistics:
# here i want to find average ages of men and women
mean(Age[female == 1])## [1] 41.87525sd(Age[female == 1])## [1] 23.87617mean(Age[!female])## [1] 39.18672sd(Age[!female])## [1] 22.90436Later you might encounter cases where you want more complicated dummy variables and want to use logical relations “and” “or” “not” (the symbols “&”, “|”, “!”) or the “>=” or multiplication or division.
As you’re going along, if you’re copying-and-pasting then you might not have had trouble, but if you’re typing then you’ve probably realized that R is persnickety - the variable Age is not equivalent to AGE nor age nor aGe… You might be wondering if there’s an easier way than copy-paste; there is. If you highlight text in the Source pane and then hit CTRL-Enter, that will run the block of code. Alt if you click the “Knit HTML” button then that will recreate the webpage as well as run all of the commands.
More Details
At the beginning you’ll be just copying commands that I give you but as you start making changes, you have to understand what’s going on.
For better or worse there is rarely just a single way of doing things, when people talk about the R language it does indeed have many features of a language. Just like a spoken language has many ways of saying “hi” so too there are lots of different ways to ask R to give means of different groups. I am showing you a few particular ways but if you look around online you’ll find others.
Remember to save your work. The computers in the lab might wipe the memory clean when you log off so back up your files. Either online (email it to yourself or upload to Google Drive or iCloud or Dropbox or whatever) and/or use a USB drive.
But be careful - you have the option of saving your code and/or your workspace; it is usually better to save the code. For example, your code might load some data (at this stage, usually data that I’ve provided). If you make changes to this data during your session, you might not want to restart with those changes. If you just save your code and re-run the program, you will start fresh. If instead you save the workspace, then you end up saving all of your scraps and changes - which eventually start to build up. So only save your workspace if you’ve run a large chunk of code that you don’t want to have to later re-run. (This is even more true if you “attach” data frames without later “detach”-ing them!)
There are 2 main file types: R script (lines of instructions to R on how to calculate some statistics) and R Markdown (which combines the instructions to R along with text around it, which is how I created this file). There are more types of course but those are the ones for now. To create files like the one you’re reading now, use “File New File R Markdown…” Give it a title and save the file with some name that you’ll remember and the extension, “.Rmd”. Keep everything in the same subdirectory or folder.
There is a bigger philosophy behind this: the idea of “reproducible research,” that if you share your .Rmd file then anybody else can run the exact same program and get the exact same results. It’s a way of convincing an audience that there was no skulduggery in how the data relates to the conclusion. It’s useful in class since it’s an easy way to submit homework or share work with your study group.
Variable Coding
Some of the PUMS variables here have a natural interpretation, for instance Age is measured in years. Actually even this has a bit of a twist, look at the histogram.
hist(Age[(Age > 45)])
A simple histogram doesn’t show it well, but a fancier version (a nonparametric kernel density), will definitely show issues:
plot(density(Age[runif(NN_obs) < 0.1], bw = "sj"))
Note, for you geeks, that if you tried this with the whole sample then the 200,000 observations would slow the computer too much, so I’ve randomly selected 10% (that’s the runif(NN_obs)<0.1 portion).
The kernel density is like a histogram with much smaller bins; in this case it shows a bit of weirdness particularly in the right, where it looks like there are suddenly a bunch of people who are 94 but nobody in between. This is due to a coding choice by the Census, where really old people are just labeled as “94” (top-coding) so it actually should be interpreted as meaning “90 or older”. So if you were to get finicky (and every good statistician is!) you might go back to the calculations of averages previously and modify them all like this, to select just those who are female and who are coded as having age less than 90. And recall that topcoding wouldn’t change the median values calculated before, which is a point in favor of that statistic
mean(Age[ (female == 1) & (Age<90) ]) ## [1] 41.13141You go make those other changes, figure out how top-coding changes the calculations of average age by gender – I’ll wait right here…
Variable Coding Again
So we were saying that some variables, like Age - ahem! – have a natural interpretation as a number.
Others are logical variables (called dummies) like female, Hispanic, or married - there is a yes/no answer that is coded 1/0. Note that if you’re creating these on your own it’s good to give names that have that sort of yes/no answer, so a variable named ‘female’ is better than one named ‘gender’ where you’d have to remember who are true and who are false.
Many variables, like PUMA, have no natural explanation at all. Here are the first codes,
str(PUMA)## num [1:196314] 3106 3106 100 4103 4103 ...You have to go to the codebook (or, in this case, my file pums_initial_recoding.r from the zip file) to find out that PUMA is “Public Use Microdata Area” where 3801 codes for Washington Heights/Inwood, 3802 is Hamilton Heights/Manhattanville/West Harlem, etc. The program will happily calculate the average value for PUMA (type in mean(PUMA) and see for yourself!) but this is a meaningless value – wtf is the average neighborhood code value!? If you want to select just people living in a particular neighborhood then you’d have to look at the list below.
| PUMA | Neighborhood |
|---|---|
| 3701 | NYC-Bronx CD 8–Riverdale, Fieldston & Kingsbridge |
| 3702 | NYC-Bronx CD 12–Wakefield, Williamsbridge & Woodlawn |
| 3703 | NYC-Bronx CD 10–Co-op City, Pelham Bay & Schuylerville |
| 3704 | NYC-Bronx CD 11–Pelham Parkway, Morris Park & Laconia |
| 3705 | NYC-Bronx CD 3 & 6–Belmont, Crotona Park East & East Tremont |
| 3706 | NYC-Bronx CD 7–Bedford Park, Fordham North & Norwood |
| 3707 | NYC-Bronx CD 5–Morris Heights, Fordham South & Mount Hope |
| 3708 | NYC-Bronx CD 4–Concourse, Highbridge & Mount Eden |
| 3709 | NYC-Bronx CD 9–Castle Hill, Clason Point & Parkchester |
| 3710 | NYC-Bronx CD 1 & 2–Hunts Point, Longwood & Melrose |
| 3801 | NYC-Manhattan CD 12–Washington Heights, Inwood & Marble Hill |
| 3802 | NYC-Manhattan CD 9–Hamilton Heights, Manhattanville & West Harlem |
| 3803 | NYC-Manhattan CD 10–Central Harlem |
| 3804 | NYC-Manhattan CD 11–East Harlem |
| 3805 | NYC-Manhattan CD 8–Upper East Side |
| 3806 | NYC-Manhattan CD 7–Upper West Side & West Side |
| 3807 | NYC-Manhattan CD 4 & 5–Chelsea, Clinton & Midtown Business District |
| 3808 | NYC-Manhattan CD 6–Murray Hill, Gramercy & Stuyvesant Town |
| 3809 | NYC-Manhattan CD 3–Chinatown & Lower East Side |
| 3810 | NYC-Manhattan CD 1 & 2–Battery Park City, Greenwich Village & Soho |
| 3901 | NYC-Staten Island CD 3–Tottenville, Great Kills & Annadale |
| 3902 | NYC-Staten Island CD 2–New Springville & South Beach |
| 3903 | NYC-Staten Island CD 1–Port Richmond, Stapleton & Mariner’s Harbor |
| 4001 | NYC-Brooklyn CD 1–Greenpoint & Williamsburg |
| 4002 | NYC-Brooklyn CD 4âBushwick |
| 4003 | NYC-Brooklyn CD 3–Bedford-Stuyvesant |
| 4004 | NYC-Brooklyn CD 2–Brooklyn Heights & Fort Greene |
| 4005 | NYC-Brooklyn CD 6–Park Slope, Carroll Gardens & Red Hook |
| 4006 | NYC-Brooklyn CD 8–Crown Heights North & Prospect Heights |
| 4007 | NYC-Brooklyn CD 16–Brownsville & Ocean Hill |
| 4008 | NYC-Brooklyn CD 5–East New York & Starrett City |
| 4009 | NYC-Brooklyn CD 18–Canarsie & Flatlands |
| 4010 | NYC-Brooklyn CD 17–East Flatbush, Farragut & Rugby |
| 4011 | NYC-Brooklyn CD 9–Crown Heights South, Prospect Lefferts & Wingate |
| 4012 | NYC-Brooklyn CD 7–Sunset Park & Windsor Terrace |
| 4013 | NYC-Brooklyn CD 10–Bay Ridge & Dyker Heights |
| 4014 | NYC-Brooklyn CD 12–Borough Park, Kensington & Ocean Parkway |
| 4015 | NYC-Brooklyn CD 14–Flatbush & Midwood |
| 4016 | NYC-Brooklyn CD 15–Sheepshead Bay, Gerritsen Beach & Homecrest |
| 4017 | NYC-Brooklyn CD 11–Bensonhurst & Bath Beach |
| 4018 | NYC-Brooklyn CD 13–Brighton Beach & Coney Island |
| 4101 | NYC-Queens CD 1–Astoria & Long Island City |
| 4102 | NYC-Queens CD 3–Jackson Heights & North Corona |
| 4103 | NYC-Queens CD 7–Flushing, Murray Hill & Whitestone |
| 4104 | NYC-Queens CD 11–Bayside, Douglaston & Little Neck |
| 4105 | NYC-Queens CD 13–Queens Village, Cambria Heights & Rosedale |
| 4106 | NYC-Queens CD 8–Briarwood, Fresh Meadows & Hillcrest |
| 4107 | NYC-Queens CD 4–Elmhurst & South Corona |
| 4108 | NYC-Queens CD 6–Forest Hills & Rego Park |
| 4109 | NYC-Queens CD 2–Sunnyside & Woodside |
| 4110 | NYC-Queens CD 5–Ridgewood, Glendale & Middle Village |
| 4111 | NYC-Queens CD 9–Richmond Hill & Woodhaven |
| 4112 | NYC-Queens CD 12–Jamaica, Hollis & St. Albans |
| 4113 | NYC-Queens CD 10–Howard Beach & Ozone Park |
| 4114 | NYC-Queens CD 14–Far Rockaway, Breezy Point & Broad Channel |
Now you’re probably thinking, isn’t there some easier way? Yes there is. R has variables called “factors” that join together the long list of codes with a separate file telling what those codes mean. Later when we do further statistics, R will know how to appropriately treat these factors. (Also it will then give an error if you calculate mean(PUMA), which is proper.)
PUMA <- as.factor(PUMA)
female <- as.factor(female)I will leave you to worry over the recoding of the other variables, because it’s good for the soul. I will show you 2 ways – the quick and dirty way, and the fancy correct way.
First the quick and dirty way.
print(levels(female))## [1] "0" "1"levels(female) <- c("male","female")Well, ways,
educ_indx <- factor((educ_nohs + 2*educ_hs + 3*educ_smcoll + 4*educ_as + 5*educ_bach + 6*educ_adv), levels=c(1,2,3,4,5,6),labels = c("No HS","HS","SmColl","AS","Bach","Adv"))These just type in the levels. But for things like PUMA, it could be a long list and might not even match every one. To do it better, we need help from a R package.
Detour on Packages
But first a bit of a detour, to mention how to use packages. R depends crucially on “packages” - that’s the whole reason that the open-source works. Some statistician invents a cool new technique, then writes up the code in R and makes it available. If you used a commercial program you’d have to wait a decade for them to update it; in R it’s here now. Also if somebody hacks a nicer or easier way to do stuff, they write it up.
So enter this into the Console,
install.packages("plyr")then
library(plyr)
levels_n <- read.csv("PUMA_levels.csv")
levels_orig <- levels(PUMA)
levels_new <- join(data.frame(levels_orig),data.frame(levels_n))## Joining by: levels_origlevels(PUMA) <- levels_new$New_LevelAlt, from R-Studio, click “Tools” then “Install Packages…” and tell it to install the package, “plyr”. That is nice if you want to see some of the packages or if you don’t quite remember the name. Then the next piece of code, library, tells the program that you want to use commands from this package.
Those commands read in a little csv file that I had made, with all of the PUMA codes, then matches the old codes with the new complete text.
Back from Detour
R will do the summary differently when it knows the variable is a factor,
summary(female)## male female
## 94400 101914summary(PUMA)## NYC-Bronx CD 8--Riverdale, Fieldston & Kingsbridge
## 915
## NYC-Bronx CD 12--Wakefield, Williamsbridge & Woodlawn
## 1251
## NYC-Bronx CD 10--Co-op City, Pelham Bay & Schuylerville
## 780
## NYC-Bronx CD 11--Pelham Parkway, Morris Park & Laconia
## 1230
## NYC-Bronx CD 3 & 6--Belmont, Crotona Park East & East Tremont
## 1358
## NYC-Bronx CD 7--Bedford Park, Fordham North & Norwood
## 804
## NYC-Bronx CD 5--Morris Heights, Fordham South & Mount Hope
## 1047
## NYC-Bronx CD 4--Concourse, Highbridge & Mount Eden
## 908
## NYC-Bronx CD 9--Castle Hill, Clason Point & Parkchester
## 1346
## NYC-Bronx CD 1 & 2--Hunts Point, Longwood & Melrose
## 1320
## NYC-Manhattan CD 12--Washington Heights, Inwood & Marble Hill
## 1133
## NYC-Manhattan CD 9--Hamilton Heights, Manhattanville & West Harlem
## 813
## NYC-Manhattan CD 10--Central Harlem
## 882
## NYC-Manhattan CD 11--East Harlem
## 839
## NYC-Manhattan CD 8--Upper East Side
## 1074
## NYC-Manhattan CD 7--Upper West Side & West Side
## 872
## NYC-Manhattan CD 4 & 5--Chelsea, Clinton & Midtown Business District
## 912
## NYC-Manhattan CD 6--Murray Hill, Gramercy & Stuyvesant Town
## 683
## NYC-Manhattan CD 3--Chinatown & Lower East Side
## 1048
## NYC-Manhattan CD 1 & 2--Battery Park City, Greenwich Village & Soho
## 960
## NYC-Staten Island CD 3--Tottenville, Great Kills & Annadale
## 1287
## NYC-Staten Island CD 2--New Springville & South Beach
## 1146
## NYC-Staten Island CD 1--Port Richmond, Stapleton & Mariner's Harbor
## 1564
## NYC-Brooklyn CD 1--Greenpoint & Williamsburg
## 1390
## NYC-Brooklyn CD 4--Bushwick
## 1073
## NYC-Brooklyn CD 3--Bedford-Stuyvesant
## 1181
## NYC-Brooklyn CD 2--Brooklyn Heights & Fort Greene
## 1284
## NYC-Brooklyn CD 6--Park Slope, Carroll Gardens & Red Hook
## 1053
## NYC-Brooklyn CD 8--Crown Heights North & Prospect Heights
## 989
## NYC-Brooklyn CD 16--Brownsville & Ocean Hill
## 1041
## NYC-Brooklyn CD 5--East New York & Starrett City
## 1439
## NYC-Brooklyn CD 18--Canarsie & Flatlands
## 2419
## NYC-Brooklyn CD 17--East Flatbush, Farragut & Rugby
## 1335
## NYC-Brooklyn CD 9--Crown Heights South, Prospect Lefferts & Wingate
## 849
## NYC-Brooklyn CD 7--Sunset Park & Windsor Terrace
## 1384
## NYC-Brooklyn CD 10--Bay Ridge & Dyker Heights
## 1473
## NYC-Brooklyn CD 12--Borough Park, Kensington & Ocean Parkway
## 1847
## NYC-Brooklyn CD 14--Flatbush & Midwood
## 1400
## NYC-Brooklyn CD 15--Sheepshead Bay, Gerritsen Beach & Homecrest
## 1625
## NYC-Brooklyn CD 11--Bensonhurst & Bath Beach
## 1982
## NYC-Brooklyn CD 13--Brighton Beach & Coney Island
## 850
## NYC-Queens CD 1--Astoria & Long Island City
## 1712
## NYC-Queens CD 3--Jackson Heights & North Corona
## 1321
## NYC-Queens CD 7--Flushing, Murray Hill & Whitestone
## 2162
## NYC-Queens CD 11--Bayside, Douglaston & Little Neck
## 1315
## NYC-Queens CD 13--Queens Village, Cambria Heights & Rosedale
## 2132
## NYC-Queens CD 8--Briarwood, Fresh Meadows & Hillcrest
## 1293
## NYC-Queens CD 4--Elmhurst & South Corona
## 955
## NYC-Queens CD 6--Forest Hills & Rego Park
## 973
## NYC-Queens CD 2--Sunnyside & Woodside
## 1101
## NYC-Queens CD 5--Ridgewood, Glendale & Middle Village
## 1821
## NYC-Queens CD 9--Richmond Hill & Woodhaven
## 1803
## NYC-Queens CD 12--Jamaica, Hollis & St. Albans
## 2427
## NYC-Queens CD 10--Howard Beach & Ozone Park
## 1344
## NYC-Queens CD 14--Far Rockaway, Breezy Point & Broad Channel
## 966
## NA's
## 126203summary(educ_indx)## No HS HS SmColl AS Bach Adv
## 62156 43109 30930 13038 26968 20113To find mean and standard deviation by neighborhood, you could use something like this,
ddply(dat_pums_NY, .(PUMA), summarize, mean = round(mean(Age), 2), sd = round(sd(Age), 2))## PUMA mean sd
## 1 100 40.55 23.95
## 2 200 42.77 23.21
## 3 300 42.71 23.13
## 4 401 41.54 23.67
## 5 402 42.63 24.53
## 6 403 42.17 23.71
## 7 500 39.32 23.30
## 8 600 40.62 23.45
## 9 701 37.00 23.40
## 10 702 40.23 22.60
## 11 703 44.47 23.67
## 12 704 43.62 23.39
## 13 800 42.03 23.66
## 14 901 43.29 24.46
## 15 902 36.54 22.64
## 16 903 35.69 22.53
## 17 904 44.86 24.06
## 18 905 41.34 23.07
## 19 906 41.02 23.93
## 20 1000 42.58 23.02
## 21 1101 42.90 23.62
## 22 1102 42.48 23.48
## 23 1201 43.26 23.14
## 24 1202 41.09 24.60
## 25 1203 42.27 23.20
## 26 1204 44.42 23.81
## 27 1205 38.78 23.95
## 28 1206 37.67 22.35
## 29 1207 42.84 23.36
## 30 1300 41.62 22.93
## 31 1400 42.27 24.46
## 32 1500 40.68 23.39
## 33 1600 43.47 24.79
## 34 1700 41.97 23.73
## 35 1801 40.84 22.36
## 36 1802 41.45 23.19
## 37 1900 41.20 22.65
## 38 2001 39.51 23.19
## 39 2002 42.19 23.50
## 40 2100 45.62 23.89
## 41 2201 41.26 24.11
## 42 2202 43.50 23.29
## 43 2203 44.76 23.90
## 44 2300 35.66 21.84
## 45 2401 43.27 24.18
## 46 2402 43.47 23.49
## 47 2500 40.64 23.90
## 48 2600 42.86 23.88
## 49 2701 43.25 23.15
## 50 2702 41.83 23.07
## 51 2801 41.38 23.08
## 52 2802 40.01 23.22
## 53 2901 40.07 23.48
## 54 2902 40.37 23.01
## 55 2903 35.10 23.54
## 56 3001 41.72 23.46
## 57 3002 41.43 24.10
## 58 3003 30.50 23.62
## 59 3101 42.62 23.13
## 60 3102 40.34 23.79
## 61 3103 41.53 24.05
## 62 3104 38.20 24.19
## 63 3105 42.67 24.42
## 64 3106 41.20 23.97
## 65 3107 41.17 24.04
## 66 3201 42.61 25.23
## 67 3202 43.93 24.26
## 68 3203 42.75 23.44
## 69 3204 42.18 23.72
## 70 3205 42.30 23.95
## 71 3206 38.88 22.93
## 72 3207 43.33 23.17
## 73 3208 42.13 23.88
## 74 3209 42.17 23.14
## 75 3210 43.27 23.62
## 76 3211 41.91 23.63
## 77 3212 44.48 24.68
## 78 3301 40.79 23.67
## 79 3302 43.54 23.73
## 80 3303 42.63 24.01
## 81 3304 39.64 23.16
## 82 3305 47.21 24.05
## 83 3306 40.03 23.10
## 84 3307 40.95 23.46
## 85 3308 39.40 22.60
## 86 3309 41.91 23.47
## 87 3310 36.43 22.30
## 88 3311 41.30 23.97
## 89 3312 42.70 23.79
## 90 3313 39.83 23.89
## 91 3701 43.95 24.62
## 92 3702 37.54 23.02
## 93 3703 43.98 24.55
## 94 3704 40.02 23.51
## 95 3705 32.33 21.78
## 96 3706 35.52 22.54
## 97 3707 33.56 21.53
## 98 3708 34.17 21.92
## 99 3709 37.12 23.35
## 100 3710 33.19 21.30
## 101 3801 41.55 22.54
## 102 3802 34.82 21.01
## 103 3803 37.39 21.14
## 104 3804 38.27 22.19
## 105 3805 43.68 23.42
## 106 3806 44.59 23.39
## 107 3807 39.93 19.40
## 108 3808 41.10 20.83
## 109 3809 40.42 21.90
## 110 3810 37.48 20.07
## 111 3901 41.78 23.31
## 112 3902 42.60 24.31
## 113 3903 38.90 23.35
## 114 4001 33.22 20.73
## 115 4002 33.79 20.74
## 116 4003 33.95 21.87
## 117 4004 37.61 20.69
## 118 4005 37.34 20.79
## 119 4006 38.39 21.91
## 120 4007 34.20 22.56
## 121 4008 35.74 22.37
## 122 4009 39.06 22.87
## 123 4010 40.11 23.41
## 124 4011 39.05 23.90
## 125 4012 35.33 21.28
## 126 4013 41.22 23.82
## 127 4014 33.95 24.41
## 128 4015 39.55 23.84
## 129 4016 41.79 24.36
## 130 4017 41.23 23.46
## 131 4018 46.05 24.51
## 132 4101 39.62 20.60
## 133 4102 37.86 22.47
## 134 4103 43.54 23.85
## 135 4104 45.08 23.46
## 136 4105 41.86 23.53
## 137 4106 40.34 23.14
## 138 4107 38.24 21.92
## 139 4108 44.41 23.32
## 140 4109 39.38 21.51
## 141 4110 40.86 22.91
## 142 4111 38.23 21.78
## 143 4112 39.17 23.33
## 144 4113 37.64 22.01
## 145 4114 38.85 24.77Although tapply would also work fine.
Here’s the 90th and 10th percentiles of wages by neighborhood,
dat_use1 <- subset(dat_pums_NY,((income_wagesal > 0) & in_NYC))ddply(dat_use1, .(PUMA), summarize, inc90 = quantile(income_wagesal,probs = 0.9), inc10 = quantile(income_wagesal,probs = 0.1), n_obs = length(income_wagesal))## PUMA inc90 inc10 n_obs
## 1 3701 122000 4800 411
## 2 3702 75000 3800 541
## 3 3703 95000 4040 353
## 4 3704 75000 5320 533
## 5 3705 53000 2700 449
## 6 3706 66600 5540 325
## 7 3707 54000 4710 378
## 8 3708 55000 6000 369
## 9 3709 69000 5000 526
## 10 3710 49800 2840 423
## 11 3801 83000 3400 563
## 12 3802 105300 3000 368
## 13 3803 100000 4000 388
## 14 3804 90000 4500 341
## 15 3805 504000 12000 562
## 16 3806 250300 10560 468
## 17 3807 282500 13300 546
## 18 3808 200000 22160 409
## 19 3809 120000 3160 489
## 20 3810 260000 8850 578
## 21 3901 100600 6200 605
## 22 3902 105000 5000 531
## 23 3903 93000 5050 696
## 24 4001 105400 9000 687
## 25 4002 60000 7000 471
## 26 4003 73000 4760 487
## 27 4004 147600 7500 715
## 28 4005 160000 11400 601
## 29 4006 90000 6130 492
## 30 4007 63000 5360 409
## 31 4008 61200 6180 539
## 32 4009 91000 6000 1124
## 33 4010 70000 7260 607
## 34 4011 68000 3900 359
## 35 4012 76000 5000 651
## 36 4013 100000 6300 681
## 37 4014 75800 4500 673
## 38 4015 80600 5880 555
## 39 4016 98200 6000 690
## 40 4017 90000 7500 870
## 41 4018 80000 7500 351
## 42 4101 85000 8480 959
## 43 4102 70000 5860 634
## 44 4103 85000 6000 970
## 45 4104 104200 5680 639
## 46 4105 85000 5000 1003
## 47 4106 90000 3960 620
## 48 4107 68200 8690 430
## 49 4108 120000 9000 492
## 50 4109 92000 7050 556
## 51 4110 85000 7300 851
## 52 4111 77000 6000 874
## 53 4112 72000 5000 1127
## 54 4113 75000 5000 610
## 55 4114 100000 4980 370You could also use table (or crosstabs) for factors with fewer items,
table(educ_indx,female)## female
## educ_indx male female
## No HS 31376 30780
## HS 20783 22326
## SmColl 15063 15867
## AS 5471 7567
## Bach 12812 14156
## Adv 8895 11218xtabs(~educ_indx + female)## female
## educ_indx male female
## No HS 31376 30780
## HS 20783 22326
## SmColl 15063 15867
## AS 5471 7567
## Bach 12812 14156
## Adv 8895 11218Want proportions instead of counts?
prop.table(table(educ_indx,female))## female
## educ_indx male female
## No HS 0.15982559 0.15678963
## HS 0.10586611 0.11372597
## SmColl 0.07672912 0.08082460
## AS 0.02786862 0.03854539
## Bach 0.06526279 0.07210897
## Adv 0.04531006 0.05714315Remember prop.table later when we do marginals.
Try it and see what happens if you use table with PUMA…
In general, R is very flexible so there are often many different ways to get the same answer. There are some people who love to debate which is best. For now just worry about learning at least one way. Later on you can go back and refine your techniques.
Sometimes attaching a dataset makes things easier. But as you get more advanced you might find it better to include the dataset name inside the function. There are advantages and disadvantages each way and some of the intro texts suggest one or the other.
If you do a lot of analysis on a particular subgroup, it might be worthwhile to create a subset of that group, so that you don’t have to always add on logical conditions. These two sets of expressions, looking at “prime-age” people, get the same results:
mean(educ_nohs[(Age >= 25)&(Age <= 55)])## [1] 0.1176959mean(educ_hs[(Age >= 25)&(Age <= 55)])## [1] 0.2455135mean(educ_smcoll[(Age >= 25)&(Age <= 55)])## [1] 0.1727465mean(educ_as[(Age >= 25)&(Age <= 55)])## [1] 0.1012049mean(educ_bach[(Age >= 25)&(Age <= 55)])## [1] 0.2126975mean(educ_adv[(Age >= 25)&(Age <= 55)])## [1] 0.1501417# alternatively
restrict1 <- as.logical((Age >= 25)&(Age <= 55))
dat_age_primeage <- subset(dat_pums_NY, restrict1)
detach()
attach(dat_age_primeage)## The following objects are masked _by_ .GlobalEnv:
##
## female, PUMA## The following objects are masked from dat_pums_NY:
##
## AfAm, Age, Amindian, ANCESTR1D, ANCESTR2D, Asian,
## below_150poverty, below_200poverty, below_povertyline,
## BUILTYR2, Commute_bus, Commute_car, Commute_other,
## Commute_rail, Commute_subway, divwidsep, educ_adv, educ_as,
## educ_bach, educ_hs, educ_nohs, educ_smcoll, female,
## foodstamps, has_AnyHealthIns, has_PvtHealthIns, HH_income,
## Hisp_Cuban, Hisp_DomR, Hisp_Mex, Hisp_PR, Hispanic, immig,
## in_Bronx, in_Brooklyn, in_Manhattan, in_Nassau, in_NYC,
## in_Queens, in_StatenI, in_Westchester, income_total,
## income_wagesal, ind_dum, Married, occ_dum, owner_cost, PERNUM,
## PUMA, race_oth, rent_cost, ROOMS, UNITSSTR, unmarried,
## veteran, white, work_fullyrmean(educ_nohs)## [1] 0.1176959mean(educ_hs)## [1] 0.2455135mean(educ_smcoll)## [1] 0.1727465mean(educ_as)## [1] 0.1012049mean(educ_bach)## [1] 0.2126975mean(educ_adv)## [1] 0.1501417detach()So you detach the original data frame and instead attach the restricted version. Then any subsequent analysis would be just done on that subset. Just remember that you’ve done this (again, this is a good reason to save the commands in a program so you can look back) otherwise you’ll wonder why you suddenly don’t have any kids in the sample!
Why All These Details?
You might be tired and bored by these details, but note that there are actually important choices to be made here, even in simply defining variables. Take the fraught American category of “race”. This data has a variable, RACED, showing how people chose to classify themselves, as ‘White,’ ‘Black,’ ‘American Indian or Alaska Native,’ Asian, various combinations, and many more codes.
Suppose you wanted to find out how many Asians are in a particular population. You could count how many people identify themselves as Asian only; you could count how many people identify as Asian in any combination. Sometimes the choice is irrelevant; sometimes it can skew the final results (e.g. the question in some areas, are there more African-Americans or more Hispanics?).
Again, there’s no “right” way to do it because there’s no science in this peculiar-but-popular concept of “race”. People’s conceptions of themselves are fuzzy and complicated; these measures are approximations.
Basics of government race/ethnicity classification
The US government asks questions about people’s race and ethnicity. These categories are social constructs, which is a fancy way of pointing out that they are based on people’s own views of themselves (influenced by how we think that other people think of us…). Currently the standard classification asks people separately about their “race” and “ethnicity” where people can pick labels from each category in any combination.
The “race” categories include: “White alone,” “Black or African-American alone,” “American Indian alone,” “Alaska Native alone,” “American Indian and Alaska Native tribes specified; or American Indian or Alaska native, not specified and no other race,” “Asian alone,” “Native Hawaiian and other Pacific Islander alone,” “Some other race alone,” or “Two or more major race groups.” Then the supplemental race categories offer more detail.
These are a peculiar combination of very general (well over 40% of the world’s population is “Asian”) and very specific (“Alaska Native alone”) representing a peculiar history of popular attitudes in the US. Only in the 2000 Census did they start to classify people in mixed races. If you were to go back to historical US Censuses from more than a century ago, you would find that the category “race” included separate entries for Irish and French and various other nationalities. Stephen J Gould has a fascinating book, The Mismeasure of Man, discussing how early scientific classifications of humans tried to “prove” which nationalities/races/groups were the smartest. Ta-Nehisi Coates notes, “racism invented race in America.”
Note that “Hispanic” is not “race” but rather ethnicity (includes various other labels such as Spanish, Latino, etc.). So a respondent could choose “Hispanic” and any race category – some choose “White,” some choose “Black,” some might be combined with any other of those complicated racial categories.
If you wanted to create a variable for those who report themselves as African-American and Hispanic, you’d use the expression (AfAm == 1) & (Hispanic == 1); sometimes stats report for non-Hispanic whites so (white == 1) & (Hispanic != 1). You can create your own classifications depending on what questions you’re investigating.
The Census Bureau gives more information here.
All of these racial categories make some people uneasy: is the government encouraging racism by recognizing these classifications? Some other governments choose not to collect race data. But that doesn’t mean that there are no differences, only that the government doesn’t choose to measure any of these differences. In the US, government agencies such as the Census and BLS don’t generally collect data on religion.
Re-Coding complicated variables from initial data
If we want more combinations of variables then we create those. Usually a statistical analysis spends a lot of time doing this sort of housekeeping - dull but necessary. It has a variety of names: data carpentry, data munging…
Educational attainment is also classified with complicated codes in this data: the original data has code 63 to mean high school diploma, 64 for a GED, 65 for less than a year of college, etc. I have transformed them into a series of dummy variables, zero/one variables for whether a person has no high school diploma, just a high school diploma, some college (but no degree), an associate’s degree, a bachelor’s degree, or an advanced degree. As with race, there is no rule that you must always do it thus. Your own analysis might be different. (Is it valid to lump together everybody who lacks a high school diploma, no matter whether they completed just first grade or up to 12th?)
That’s the whole point of learning to do the data work for yourself: you can see all of the little decisions that go into creating a conclusion. Some conclusions might be fragile so a tiny decision about coding could change everything; other conclusions are robust to deviations. You must find out.
De-bugging
Without a doubt, programming is tough. In R or with any other program, it is frustrating and complicated and difficult to do it the first few times. Some days it feels like a continuous battle just to do the simplest thing! Keep going despite that, keep working on it.
Your study group will be very helpful of course.
Often a google search of the error message helps. If you’ve isolated the error and read the help documentation on that command, then you’re on your way to solving the problem on your own.
If you have troubles that you can’t solve, email me for help. But try to narrow down your question: if you run 20 lines of code that produce an error, is there a way to reproduce the error after just 5 lines? What if you did the same command on much simpler data, would it still cause an error? Sending emails like “I have a problem with errors” might be cathartic but is not actually useful to anyone. If you email me with the minimal code that recreates the error, along with the text of the error and/or a screenshot, then that will help more.
Do it
The first homework assignment asks you to start working on these questions. Begin by running the code that I give here, just to see if you can replicate my results. Then start asking more questions about the data - there is so much info there! Immigration/ancestry status, how they commute, health insurance, poverty, income, owner or renter (and how much they pay!), all sorts of info. Have some fun.