-
Notifications
You must be signed in to change notification settings - Fork 72
/
cm101.Rmd
796 lines (600 loc) · 33.5 KB
/
cm101.Rmd
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
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
# (1) Review of STAT 545
```{r include=FALSE}
knitr::opts_chunk$set(echo=TRUE, warning=FALSE, message=FALSE)
```
## Today's Agenda
- Part 1: Introductions and course overview (10 mins)
- Two parts to the course: reproducible workflows and dashboards
- Course grading elements
- Group projects
- Optional "lab"
- Office hours
- What to expect in STAT 547
- Part 2: Group projects in STAT 547 (10 mins)
- Partners
- Choosing a research question
- What's due this week
- Part 3: Review of STAT545 (45 mins)
- tidy data
- ggplot
- dplyr
- File I/O
- Other stuff
## Part 1: Introductions and course overview (15 mins)
STAT 547 can roughly be divided into two main ideas or themes: **reproducible workflows** and **interactive dashboards**.
In the first part of the course we will start working with writing functions in R scripts and think a little bit about functional programming.
Then with some new tools at our disposal, we will create "analysis pipelines" to automate and streamline the process going from loading raw data, all the way to creating a final report with the results of a regression analysis.
Next we'll switch gears and use our functional programming skills to create a dashboard app with multiple interactive components and multi-page/tab layouts.
Finally we will end with deploying our dashboard to a free cloud service so that it is publicly available!
It is my hope that you will use the things we learn in this course to inform your own data analysis pipelines in your PhD and Masters research projects.
To that end, if you ever have questions about these tools could be used in your projects, feel free to drop by my office hours and I'd be happy to discuss it!
### Important links
- All course announcements, discussions, and participation worksheets will be posted: [https://github.com/STAT547-UBC-2019-20/Discussions](https://github.com/STAT547-UBC-2019-20/Discussions)
- The STAT 547 GitHub organization can be accessed here: [https://github.com/STAT547-UBC-2019-20](https://github.com/STAT547-UBC-2019-20)
- All assignment/milestone/participation submissions are to be done on Canvas (submit URL of tagged release in your repo)
### Course structure
- (10%) Participation: Each week there will be two 90 minute classes closely following this guidebook. I hope to begin each meeting with a short introduction to provide some motivation for the topic of the day. We will then do coding exercises together. These will be due at the **end of the course* in your `participation` repository. Each file will be graded as complete (1), half-complete (0.5), or missing (0). Of the 12 class meetings, your 10 best will be used.
- (20%) Assignments: These will be short assignments designed to help you practice material from the current week. Your lowest assignment score will be dropped (so 5 of 6 assignments will count for your grade). Links to your `assignment_X` repo should be submitted on Canvas and will be due on **Saturdays at 6 PM**.
- (60%) Team projects: There will be two team projects in this course (each equally weighted), one building an analysis pipeline and the second creating a dashboard using DashR. You will stay with the same groups for the duration of the course and both partners will receive the same grade. Though the projects are separate, they will use the same dataset and should be considered as extensions of each other. Links to your tagged release in your group project rep (`group_##_SOME_RELATED_NAME`) repo should be submitted on Canvas and will be due on **Saturdays at 6 PM**.
- (10%) Teamwork: The teamwork document will be due after the end of each project and is your opportunity to reflect on your group dynamics and share with the teaching team your successess and challenges. It will constitute a 4 point rating scale on a few categories as well as an area for written feedback. A template teamwork document will be provided to you. This is not intended to be a long written exercise (unless you choose to do so, we will read it all).
I have added a 48-hour grace period for all submissions - this means that should "life happen", you will still be able to submit your assignments for full credit up to 48 hours after the due date.
You do not need to inform the teaching team if you need to use the grace period.
After the 48 hour grace period, solutions will be released and no submissions will be accepted - no exceptions.
### Office hours and optional "lab"
Each week there will be about 3 office hours scheduled.
Check the [course website](https://stat545.stat.ubc.ca/officehours/) for details on the office hour schedule.
### Auditing students
You're more than welcome to audit the course!
However, given the nature of this course, it will NOT mean any reduction in workload.
Auditing students will be:
1) Graded on a Pass/Fail basis for each assignment/course activity
2) Exempted from the peer review portion (assignment) of the course
3) Doing the team project as an individual project (with no other change in scope or requirements)
If you would like to audit the course, you *must* submit the request by 11:59 PM on Thursday March 5th, 2020. After this - without exceptions - students will not be permitted to switch to (or from) audit status.
## Part 2: Group projects in STAT 547
### Partners
Your partners for the group projects will be randomly assigned using an algorithm with no human intervention.
Teams will be posted on Canvas.
There will be a "milestone" assignment due each week where you will be expected to work with your partner collaboratively in a single GitHub repository using the workflows you started learning about in STAT 545.
The skills you learn in project courses working with git and GitHub will be essential for your future career - especially if it involves any amount of code.
### Milestone 1
In the first milestone, your task is to introduce yourself to your project partner, exchange contact information, choose a dataset, and do some preliminary EDA on the dataset.
We will talk about formulating a research question later on in the course, but keep in mind you will be doing a linear regression analysis to answer a particular research question.
Keep in mind that you will be working with the dataset you choose for the remainder of the course, so I've provided some [guidelines](/evaluation/projects) for choosing a dataset.
Best of luck - this project course will be tonnes of fun !
## Part 3: Review of STAT545 (45 mins)
## Participation worksheet
- Download the participation worksheet from [here](https://github.com/STAT547-UBC-2019-20/Discussions/blob/master/participation/cm101/cm101_participation.Rmd)
- Add it to your local clone of your participation repo (don't forget to ["accept"](https://classroom.github.com/a/uYT4-Doe) your participation repo first)
- Today's participation activities are all at the end, in the section "Do I need more practice with the tidyverse?"
### Introduction
This section is adapted from Chapter 3 of the [tidynomicon textbook by Greg Wilson](https://gvwilson.github.io/tidynomicon/).
### The Tidyverse {#tidyverse}
There is no point learning a language designed for data manipulation if you do not then bend data to your will.
This chapter therefore looks at how to do the things that R was summoned—er, designed—to do.
### Learning Objectives
- Install and load packages in R.
- Read CSV data with R.
- Explain what a tibble is and how tibbles related to data frames and matrices.
- Describe how `read_csv` infers data types for columns in tabular datasets.
- Name and use three functions for inspects tibbles.
- Select subsets of tabular data using column names, scalar indices, ranges, and logical expressions.
- Explain the difference between indexing with `[` and with `[[`.
- Name and use four functions for calculating aggregate statistics on tabular data.
- Explain how these functions treat `NA` by default, and how to change that behavior.
### How do I read data?
We begin by looking at the file [`infant_hiv.csv`](https://github.com/gvwilson/tidynomicon/raw/master/results/infant_hiv.csv), a tidied version of data on the percentage of infants born to women with HIV who received an HIV test themselves within two months of birth.
The original data comes from the UNICEF site at <https://data.unicef.org/resources/dataset/hiv-aids-statistical-tables/>,
and this file contains:
```
country,year,estimate,hi,lo
AFG,2009,NA,NA,NA
AFG,2010,NA,NA,NA
...
AFG,2017,NA,NA,NA
AGO,2009,NA,NA,NA
AGO,2010,0.03,0.04,0.02
AGO,2011,0.05,0.07,0.04
AGO,2012,0.06,0.08,0.05
...
ZWE,2016,0.71,0.88,0.62
ZWE,2017,0.65,0.81,0.57
```
The actual file has many more rows (and no ellipses).
It uses `NA` to show missing data rather than (for example) `-`, a space, or a blank,
and its values are interpreted as follows:
| Header | Datatype | Description |
|----------|-----------|---------------------------------------------|
| country | char | ISO3 country code of country reporting data |
| year | integer | year CE for which data reported |
| estimate | double/NA | estimated percentage of measurement |
| hi | double/NA | high end of range |
| lo | double/NA | low end of range |
We will import this data first by installing (if necessary, using `install.packages("tidyverse")`) the tidyverse collection of packages and then calling the `read_csv` function.
```{r library-install, eval = TRUE, message = TRUE, warning = TRUE}
library(tidyverse)
```
Loading the tidyverse gives us eight packages (ggplot2, tibble, `tidyr`, `readr`, `purr`, `dplyr`, `stringr`, `forcats`.
One of those, `dplyr`, defines two functions that mask standard functions in R with the same names.
If we need the originals, we can always get them with their fully-qualified names referring to the original package they belong to: `stats::filter` and `stats::lag`.
Firas' note: I find these errors and messages quite irritating so I have made it so they do not appear unless I explicitly ask it to (which I did in this case for illustrative purposes).
At any time, we can call `sessionInfo()` to find out what versions of which packages we have loaded, along with the version of R we're using and some other useful information:
```{r}
sessionInfo()
```
Once we have the tidyverse loaded, reading the file looks remarkably like reading the file:
```{r r-read-csv}
infant_hiv <- read_csv('https://github.com/gvwilson/tidynomicon/raw/master/results/infant_hiv.csv')
```
R's `read_csv` tells us a little about what it has done.
In particular, it guesses the data types of columns based on the first thousand values
and then tells us what types it has inferred.
(In a better universe, people would habitually use the first *two* rows of their spreadsheets for name *and units*, but we do not live there.)
We can now look at what `read_csv` has produced:
```{r show-tibble}
infant_hiv
```
This is a tibble, which is the tidyverse's enhanced version of R's `data.frame`.
It organizes data into named columns, each having one value for each row.
In statistical terms, the columns are the variables being observed and the rows are the actual observations.
### How do I inspect data?
We often have a quick look at the content of a table to remind ourselves what it contains.
R has `head` and `tail` functions:
```{r r-show-head}
head(infant_hiv)
```
```{r r-show-tail}
tail(infant_hiv)
```
Let's have a closer look at that last command's output:
```{r r-tail-tibble, paged.print=FALSE}
tail(infant_hiv)
```
Note that the row numbers printed by `tail` are relative to the output, not absolute to the table.
What about overall information?
```{r data-summary}
summary(infant_hiv)
```
Your display of R's summary may or may not wrap, depending on how large a screen the older acolytes have allowed you.
### How do I index rows and columns?
To access columns in R, you can use an attribute name:
```{r r-string-subscript}
infant_hiv['estimate']
```
However, this syntax `infant_hiv$estimate` provides all the data:
```{r r-dollar-subscript, output.lines=NA}
infant_hiv$estimate
```
Again, note that the boxed number on the left is the start index of that row.
What about single values?
We have:
```{r r-individual-element}
infant_hiv$estimate[12]
```
Ah—everything in R is a vector,
so we get a vector of one value as an output rather than a single value.
```{r r-len-individual-element}
length(infant_hiv$estimate[12])
```
And yes, ranges work:
```{r r-range-estimate}
infant_hiv$estimate[6:15]
```
Note that the upper bound inclusive in R.
Note also that nothing prevents us from selecting a range of rows that spans several countries,
which is why selecting by row number is usually a sign of innocence, insouciance, or desperation.
We can select by column number as well.
In R, a single index is interpreted as the column index:
```{r single-index-is-col}
infant_hiv[1]
```
But notice that the output is not a vector, but another tibble (i.e., a table with N rows and one column).
This means that adding another index does column-wise indexing on that tibble:
```{r double-index-of-tibble}
infant_hiv[1][1]
```
How then are we to get the first mention of Afghanistan?
The answer is to use double square brackets to strip away one level of structure:
```{r double-square-on-tibble, output.lines=NA}
infant_hiv[[1]]
```
This is now a plain old vector, so it can be indexed with single square brackets:
```{r double-square-then-single}
infant_hiv[[1]][1]
```
But that too is a vector, so it can of course be indexed as well (for some value of "of course"):
```{r double-single-single}
infant_hiv[[1]][1][1]
```
Thus, `data[1][[1]]` produces a tibble, then selects the first column vector from it, so it still gives us a vector.
*This is not madness.*
It is merely...differently sane.
> **Subsetting Data Frames**
>
> When we are working with data frames (including tibbles),
> subsetting with a single vector selects columns, not rows,
> because data frames are stored as lists of columns.
> This means that `df[1:2]` selects two columns from `df`.
> However, in `df[2:3, 1:2]`, the first index selects rows, while the second selects columns.
### How do I calculate basic statistics?
What is the average estimate?
We start by grabbing that column for convenience:
```{r r-average-index}
estimates <- infant_hiv$estimate
length(estimates)
```
```{r r-estimates-mean}
mean(estimates)
```
The void is always there, waiting for us...
Let's fix this in R first by telling `mean` to drop NAs:
```{r r-remove-na}
mean(estimates, na.rm = TRUE)
```
Many functions in R use `na.rm` to control whether `NA`s are removed or not.
(Remember, the `.` character is just another part of the name)
R's default behavior is to leave `NA`s in, and then to include them in aggregate computations.
A good use of aggregation is to check the quality of the data.
For example, we can ask if there are any records where some of the estimate, the low value, or the high value are missing,
```{r r-check-null}
any(is.na(infant_hiv$hi) != is.na(infant_hiv$lo))
```
**Aside:** In R, there's a handy package called `glue`.
The `glue` package is very useful for "gluing" strings to data in R.
Let's load it to see how it works:
```{r r-min-max-std}
library(glue)
print(glue("min {min(estimates, na.rm = TRUE)}"))
print(glue("max {max(estimates, na.rm = TRUE)}"))
print(glue("sd {sd(estimates, na.rm = TRUE)}"))
```
### How do I filter data?
By "filtering", we mean "selecting records by value".
The simplest approach is to use a vector of logical values to keep only the values corresponding to `TRUE`:
```{r r-simple-filter}
maximal <- estimates[estimates >= 0.95]
length(maximal)
```
```{r r-maximal, output.lines=NA}
maximal
```
It appears that R has kept the unknown values in order to highlight just how little we know.
More precisely, wherever there was an `NA` in the original data, there is an `NA` in the logical vector and hence, an `NA` in the final vector.
Let us then turn to `which` to get a vector of indices at which a vector contains `TRUE`.
This function does not return indices for `FALSE` or `NA`:
```{r r-which}
which(estimates >= 0.95)
```
And as a quick check:
```{r r-length-which}
length(which(estimates >= 0.95))
```
So now we can index our vector with the result of the `which`:
```{r r-maximal-which}
maximal <- estimates[which(estimates >= 0.95)]
maximal
```
But should we do this?
Those `NA`s are important information, and should not be discarded so blithely.
What we should *really* be doing is using the tools the tidyverse provides rather than clever indexing tricks.
These behave consistently across a wide scale of problems and encourage use of patterns that make it easier for others to understand our programs.
### How do I write tidy code?
The six basic data transformation operations in the tidyverse are:
- `filter`: choose observations (rows) by value(s)
- `arrange`: reorder rows
- `select`: choose variables (columns) by name
- `mutate`: derive new variables from existing ones
- `group_by`: define subsets of rows for further processing
- `summarize`: combine many values to create a single new value
`filter(tibble, ...criteria...)` keeps rows that pass all of the specified criteria:
```{r filter-as-function}
filter(infant_hiv, lo > 0.5)
```
Notice that the expression is `lo > 0.5` rather than `"lo" > 0.5`.
The latter expression would return the entire table because the string `"lo"` is greater than the number 0.5 everywhere.
But how is it that the name `lo` can be used on its own?
It is the name of a column, but there is no variable called `lo`.
The answer is that R uses lazy evaluation: function arguments aren't evaluated until they're needed, so the function `filter` actually gets the expression `lo > 0.5`, which allows it to check that there's a column called `lo` and then use it appropriately.
It may seem strange at first, but it is much tidier than `filter(data, data$lo > 0.5)` or `filter(data, "lo > 0.5")`.
We can make data anlaysis code more readable by using the pipe operator `%>%`:
```{r filter-in-pipe}
infant_hiv %>% filter(lo > 0.5)
```
This may not seem like much of an improvement,
but neither does a Unix pipe consisting of `cat filename.txt | head`.
What about this?
```{r filter-complex}
filter(infant_hiv, (estimate != 0.95) & (lo > 0.5) & (hi <= (lo + 0.1)))
```
It uses the vectorized "and" operator `&` twice,
and parsing the condition takes a human being at least a few seconds.
Its pipelined equivalent is:
```{r filter-complex-pipe}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.1))
```
Breaking the condition into stages like this often makes reading and testing much easier, and encourages incremental write-test-extend development.
Let's increase the band from 10% to 20%, break the line the way the tidyverse style guide recommends to make the operations easier to spot, and order by `lo` in descending order:
```{r filter-arrange}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.2)) %>%
arrange(desc(lo))
```
We can now `select` the three columns we care about:
```{r filter-arrange-select}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.2)) %>%
arrange(desc(lo)) %>%
select(year, lo, hi)
```
Once again, we are using the unquoted column names `year`, `lo`, and `hi` and letting R's lazy evaluation take care of the details for us.
Rather than selecting these three columns, we can select *out* the columns we're not interested in by negating their names.
This leaves the columns that are kept in their original order, rather than putting `lo` before `hi`, which won't matter if we later select by name, but *will* if we ever want to select by position:
```{r select-out}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.2)) %>%
arrange(desc(lo)) %>%
select(-country, -estimate)
```
Giddy with power, we now add a column containing the difference between the low and high values.
This can be done using either `mutate`, which adds new columns to the end of an existing tibble, or with `transmute`, which creates a new tibble containing only the columns we explicitly ask for.
(There is also a function `rename` which simply renames columns.)
Since we want to keep `hi` and `lo`, we decide to use `mutate`:
```{r mutate-new-column}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.2)) %>%
arrange(desc(lo)) %>%
select(-country, -estimate) %>%
mutate(difference = hi - lo)
```
Does the difference between high and low estimates vary by year?
To answer that question,
we use `group_by` to group records by value and then `summarize` to aggregate within groups.
We might as well get rid of the `arrange` and `select` calls in our pipeline at this point (since we're not using them) and count how many records contributed to each aggregation using `n()`:
```{r summarize-and-count}
infant_hiv %>%
filter(estimate != 0.95) %>%
filter(lo > 0.5) %>%
filter(hi <= (lo + 0.2)) %>%
mutate(difference = hi - lo) %>%
group_by(year) %>%
summarize(count = n(), ave_diff = mean(year))
```
Elegant, right?
### How do I model my data?
Tidying up data can be as calming and rewarding in the same way as knitting or rearranging the specimen jars on the shelves in your dining room-stroke-laboratory.
Eventually, though, people want to do some statistics.
The simplest tool for this in R is `lm`, which stands for "linear model".
Given a formula and a data set, it calculates coefficients to fit that formula to that data:
```{r simple-formula}
model <- lm(estimate ~ lo, data = infant_hiv)
model
```
This is telling us that `estimate` is more-or-less equal to `0.0421 + 1.0707 * lo`.
The `~` symbol is used to separate the left and right sides of the equation, and as with all things tidyverse, lazy evaluation allows us to use variable names directly.
Recall another handy R package `broom` that helps us work with the output of `lm` in the tidyverse way:
```{r}
library(broom)
tidy(model)
```
You'll note that the same information is now presented to usin a more familiar form (a tibble).
The row names have been moved into a column called term, and the column names are consistent and can be accessed using $.
In case you are interested in the fitted values and residuals for each of the original points in the regression.
For this, use `augment`:
```{r}
augment(model)
```
Note that each of the new columns begins with a '.' to avoid overwriting any of the original columns.
Finally, several summary statistics are computed for the entire regression, such as R^2 and the F-statistic.
These can be accessed with the `glance` function:
```{r}
glance(model)
```
In fact, it lets us write much more complex formulas involving functions of multiple variables.
For example, we can regress `estimate` against the square roots of `lo` and `hi` (though there is no sound statistical reason to do so):
```{r complex-formula}
model2 <- lm(estimate ~ sqrt(lo) + sqrt(hi), data = infant_hiv)
model2
```
Let's take amoment to describe how the `+` is overloaded in formulas.
The formula `estimate ~ lo + hi` does *not* mean "regress `estimate` against the sum of `lo` and `hi`",
but rather, "regress `estimate` against the two variables `lo` and `hi`":
```{r double-regression}
model3 <- lm(estimate ~ lo + hi, data = infant_hiv)
model3
```
If we want to regress `estimate` against the average of `lo` and `hi`
(i.e., regress `estimate` against a single calculated variable instead of against two variables)
we need to create a temporary column:
```{r regress-temporary}
infant_hiv %>%
mutate(ave_lo_hi = (lo + hi)/2) %>%
lm(estimate ~ ave_lo_hi, data = .)
```
Here, the call to `lm` is using the variable `.` to mean "the data coming in from the previous stage of the pipeline".
Most of the functions in the tidyverse use this convention so that data can be passed to a function that expects it in a position other than the first.
### How do I create a plot?
Human being always want to see the previously unseen, though they are not always glad to have done so.
The most popular tool for doing this in R is `ggplot2`, which implements and extends the patterns for creating charts described in @Wilk2005.
Every chart it creates has a geometry that controls how data is displayed and a mapping that controls how values are represented geometrically.
For example, these lines of code create a scatter plot showing the relationship between `lo` and `hi` values in the infant HIV data:
```{r basic-plot}
ggplot(infant_hiv) +
geom_point(mapping = aes(x = lo, y = hi)) +
theme_bw()
# Missing proper axed labels and titles
```
Looking more closely:
- The function `ggplot` creates an object to represent the chart with `infant_hiv` as the underlying data.
- `geom_point` specifies the geometry we want (points).
- Its `mapping` argument is assigned an [aesthetic](glossary.html#aesthetic)
that specifies `lo` is to be used as the `x` coordinate and `hi` is to be used as the `y` coordinate.
- The `theme_bw()` is added to the plot to remove the default grey background
- The elements of the chart are combined with `+` rather than `%>%` for historical reasons.
Let's create a slightly more appealing plot by dropping NAs, making the points semi-transparent, and colorizing them according to the value of `estimate`:
```{r plot-after-drop}
infant_hiv %>%
drop_na() %>%
ggplot(mapping = aes(x = lo, y = hi, color = estimate)) +
geom_point(alpha = 0.5) +
xlim(0.0, 1.0) + ylim(0.0, 1.0) +
theme_bw() +
labs(x = 'Low estimate',
y = 'High estimate',
title = '% of infants born to women with HIV ... test themselves w/in two months of birth')
```
We set the transparency `alpha` outside the aesthetic because its value is constant for all points.
If we set it inside `aes(...)`, we would be telling ggplot2 to set the transparency according to the value of the data.
We specify the limits to the axes manually with `xlim` and `ylim` to ensure that ggplot2 includes the upper bounds.
Without this, all of the data would be shown, but the upper label "1.00" would be omitted.
This plot immediately shows us that we have some outliers.
There are far more values with `hi` equal to 0.95 than it seems there ought to be, and there are eight points running up the left margin that seem troubling as well.
Let's create a new tibble that doesn't have these:
```{r plot-remove-outliers}
infant_hiv %>%
drop_na() %>%
filter(hi != 0.95) %>%
filter(!((lo < 0.10) & (hi > 0.25))) %>%
ggplot(mapping = aes(x = lo, y = hi, color = estimate)) +
geom_point(alpha = 0.5) +
xlim(0.0, 1.0) + ylim(0.0, 1.0) +
theme_bw() +
labs(x = 'Low estimate',
y = 'High estimate',
title = '% of infants born to women with HIV ... test themselves w/in two months of birth')
```
We can add the fitted curve by including another geometry called `geom_smooth`:
```{r plot-with-fit}
infant_hiv %>%
drop_na() %>%
filter(hi != 0.95) %>%
filter(!((lo < 0.10) & (hi > 0.25))) %>%
ggplot(mapping = aes(x = lo, y = hi)) +
geom_point(mapping = aes(color = estimate), alpha = 0.5) +
geom_smooth(method = lm, color = 'red') +
xlim(0.0, 1.0) + ylim(0.0, 1.0) +
theme_bw() +
labs(x = 'Low estimate',
y = 'High estimate',
title = '% of infants born to women with HIV ... test themselves w/in two months of birth')
```
But wait: why is this complaining about missing values?
Some online searches lead to the discovery that `geom_smooth` adds virtual points to the data for plotting purposes, some of which lie outside the range of the actual data, and that setting `xlim` and `ylim` then truncates these.
(Remember, R is differently sane...)
The safe way to control the range of the data is to add a call to `coord_cartesian`, which effectively zooms in on a region of interest:
```{r plot-cartesian}
infant_hiv %>%
drop_na() %>%
filter(hi != 0.95) %>%
filter(!((lo < 0.10) & (hi > 0.25))) %>%
ggplot(mapping = aes(x = lo, y = hi)) +
geom_point(mapping = aes(color = estimate), alpha = 0.5) +
geom_smooth(method = lm, color = 'red') +
coord_cartesian(xlim = c(0.0, 1.0), ylim = c(0.0, 1.0)) +
theme_bw() +
labs(x = 'Low estimate',
y = 'High estimate',
title = '% of infants born to women with HIV ... test themselves w/in two months of birth')
```
### Do I need more practice with the tidyverse?
Absolutely!
1. First, accept your [`participation` assignment using GitHub classroom](https://classroom.github.com/a/uYT4-Doe)
2. Navigate to the `cm101` folder. Create a new subfolder inside here called `data`.
3. In the `cm101_participation.Rmd` file, load the `here` package (recall it is used to construct paths without specifying user-specific paths).
```{r fake-load-libraries, eval=FALSE}
library(tidyverse)
library(here)
```
4. Download the [person.csv](https://github.com/gvwilson/tidynomicon/raw/master/data/person.csv) file into the `data` folder you created.
```{r download_person}
readr::write_csv(readr::read_csv("https://github.com/gvwilson/tidynomicon/raw/master/data/person.csv"),here::here("data", "person.csv"))
```
4. Next, use `here::here` to construct a path to a file and `readr::read_csv` to read the `person.csv` file:
```{r read-survey-data}
path = here::here("data", "person.csv")
person <- readr::read_csv(path)
```
We don't need to write out fully-qualified names—`here` and `read_csv` will do—but we will use them to make it easier to see what comes from where.
5. Next, have a look at the tibble `person`, which contains some basic information about a group of foolhardy scientists who ventured into the Antarctic in the 1920s and 1930s in search of things best left undisturbed:
```{r show-person}
person
```
6. How many rows and columns does this tibble contain?
```{r count-rows}
nrow(person)
```
```{r count-cols}
ncol(person)
```
(These names don't have a package prefix because they are built in.)
7. Let's show that information in a slightly nicer way using `glue` to insert values into a string and `print` to display the result:
```{r use-glue}
print(glue::glue("person has {nrow(person)} rows and {ncol(person)} columns"))
```
8. Use the function `paste` to combine the elements of a vector to display several values.`colnames` gives us the names of a tibble's columns, and `paste`'s `collapse` argument tells the function to use a single space to separate concatenated values:
```{r use-colnames-and-paste}
print(glue::glue("person columns are {paste(colnames(person), collapse = ' ')}"))
```
9. Time for some data manipulation. Let's get everyone's family and personal names:
```{r select-by-name}
dplyr::select(person, family_name, personal_name)
```
10. Then filter that list to keep only those that come in the first half of the alphabet:
```{r filter-with-two-conditions}
dplyr::select(person, family_name, personal_name) %>%
dplyr::filter(family_name < "N")
```
11. It would be more consistent to rewrite this as:
```{r filter-consistently}
person %>%
dplyr::select(family_name, personal_name) %>%
dplyr::filter(family_name < "N")
```
12. Let's add a column that records the lengths of family names:
```{r mutate-name-length}
person %>%
dplyr::mutate(name_length = stringr::str_length(family_name))
```
13. Then arrange in descending order:
```{r mutate-and-arrange}
person %>%
dplyr::mutate(name_length = stringr::str_length(family_name)) %>%
dplyr::arrange(dplyr::desc(name_length))
```
### Key Points
- `install.packages('name')` installs packages.
- `library(name)` (without quoting the name) loads a package.
- `library(tidyverse)` loads the entire collection of tidyverse libraries at once.
- `read_csv(filename)` reads CSV files that use the string 'NA' to represent missing values.
- `read_csv` infers each column's data types based on the first thousand values it reads.
- A tibble is the tidyverse's version of a data frame, which represents tabular data.
- `head(tibble)` and `tail(tibble)` inspect the first and last few rows of a tibble.
- `summary(tibble)` displays a summary of a tibble's structure and values.
- `tibble$column` selects a column from a tibble, returning a vector as a result.
- `tibble['column']` selects a column from a tibble, returning a tibble as a result.
- `tibble[,c]` selects column `c` from a tibble, returning a tibble as a result.
- `tibble[r,]` selects row `r` from a tibble, returning a tibble as a result.
- Use ranges and logical vectors as indices to select multiple rows/columns or specific rows/columns from a tibble.
- `tibble[[c]]` selects column `c` from a tibble, returning a vector as a result.
- `min(...)`, `mean(...)`, `max(...)`, and `std(...)` calculates the minimum, mean, maximum, and standard deviation of data.
- These aggregate functions include `NA`s in their calculations, and so will produce `NA` if the input data contains any.
- Use `func(data, na.rm = TRUE)` to remove `NA`s from data before calculations are done (but make sure this is statistically justified).
- `filter(tibble, condition)` selects rows from a tibble that pass a logical test on their values.
- `arrange(tibble, column)` or `arrange(desc(column))` arrange rows according to values in a column (the latter in descending order).
- `select(tibble, column, column, ...)` selects columns from a tibble.
- `select(tibble, -column)` selects *out* a column from a tibble.
- `mutate(tibble, name = expression, name = expression, ...)` adds new columns to a tibble using values from existing columns.
- `group_by(tibble, column, column, ...)` groups rows that have the same values in the specified columns.
- `summarize(tibble, name = expression, name = expression)` aggregates tibble values (by groups if the rows have been grouped).
- `tibble %>% function(arguments)` performs the same operation as `function(tibble, arguments)`.
- Use `%>%` to create pipelines in which the left side of each `%>%` becomes the first argument of the next stage.