1019hw3

data input

pacman::p_load(foreign)
dta <- read.dta("bodyfat.dta", 
                 convert.dates = TRUE,
                 convert.factors = TRUE,
                 missing.type = FALSE,
                 convert.underscore = FALSE,
                 warn.missing.labels = TRUE)
str(dta)

## 'data.frame':    2273 obs. of  4 variables:
##  $ id     : num  1 1 1 2 2 2 3 3 3 4 ...
##  $ sex    : Factor w/ 2 levels "female","male": 2 2 2 1 1 1 1 1 1 1 ...
##  $ age    : num  11 13 15 11 13 15 11 13 15 11 ...
##  $ bodyfat: num  4 6.2 10.5 8.1 10.4 ...
##  - attr(*, "datalabel")= chr ""
##  - attr(*, "time.stamp")= chr "21 May 2011 13:30"
##  - attr(*, "formats")= chr [1:4] "%9.0g" "%9.0g" "%9.0g" "%9.0g"
##  - attr(*, "types")= int [1:4] 254 254 254 254
##  - attr(*, "val.labels")= chr [1:4] "" "labsex" "" ""
##  - attr(*, "var.labels")= chr [1:4] "group(u)" "" "" ""
##  - attr(*, "version")= int 8
##  - attr(*, "label.table")=List of 1
##   ..$ labsex: Named int [1:2] 1 2
##   .. ..- attr(*, "names")= chr [1:2] "female" "male"

head (dta)

##   id    sex age bodyfat
## 1  1   male  11     4.0
## 2  1   male  13     6.2
## 3  1   male  15    10.5
## 4  2 female  11     8.1
## 5  2 female  13    10.4
## 6  2 female  15    15.2

ggplot

library(ggplot2)
ggplot(dta, aes(age, bodyfat,
                  group = id,
                  color = id)) +
  geom_point()+
  stat_summary(fun = mean, geom = "line") +
  #stat_summary(fun = mean, geom = "point") +
  #stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.3) +
  #scale_shape_manual(values = c(1, 2)) +
  facet_wrap( ~ sex)+
  labs(x = "Age (in years)", 
       y = "Body fat (in %)") +
  theme(legend.justification = c(1, 1), 
        legend.position = c(1, 1),
        legend.background = element_rect(fill = "white",color = "black"))

descriptive table

pacman::p_load(dplyr, furniture)

dta %>% 
  dplyr::group_by(sex) %>%
  dplyr::select(starts_with("bodyfat")) %>%
  furniture::table1(digits=2, total=FALSE, test=F, output="html")

## Adding missing grouping variables: `sex`

## Using dplyr::group_by() groups: sex

	female	male
	n = 1045	n = 1228
bodyfat
	10.04 (2.56)	8.84 (2.49)