Hw7
Max Bauer
2022-11-03
library(tidyr)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(OpenMx)## Warning: package 'OpenMx' was built under R version 4.2.2data(twinData)
twinData <- as_tibble(twinData)twinData## # A tibble: 3,808 × 16
## fam age zyg part wt1 wt2 ht1 ht2 htwt1 htwt2 bmi1 bmi2
## <int> <int> <int> <int> <int> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 21 1 2 58 57 1.7 1.7 20.1 19.7 21.0 20.9
## 2 2 24 1 2 54 53 1.63 1.63 20.3 19.9 21.1 21.0
## 3 3 21 1 2 55 50 1.65 1.68 20.2 17.7 21.0 20.1
## 4 4 21 1 2 66 76 1.57 1.65 26.8 27.9 23.0 23.3
## 5 5 19 1 2 50 48 1.61 1.63 19.3 18.1 20.7 20.3
## 6 6 26 1 2 60 60 1.60 1.57 23.4 24.3 22.1 22.3
## 7 7 23 1 2 65 65 1.75 1.77 21.2 20.7 21.4 21.2
## 8 8 29 1 2 40 39 1.56 1.53 16.4 16.7 19.6 19.7
## 9 9 24 1 2 60 57 1.76 1.77 19.4 18.2 20.7 20.3
## 10 10 28 1 2 76 64 1.7 1.73 26.3 21.4 22.9 21.4
## # … with 3,798 more rows, and 4 more variables: cohort <chr>, zygosity <fct>,
## # age1 <int>, age2 <int>library(ggplot2)ggplot(twinData, aes(x=ht1, y=ht2)) + geom_point()## Warning: Removed 141 rows containing missing values (geom_point).
Yes because there is a positive coorelation between the two.
twinData2 <- ggplot(twinData, aes(x=ht1, y=ht2)) + geom_point()
twinData2 + facet_wrap(~cohort)## Warning: Removed 141 rows containing missing values (geom_point).
twinData3 <- twinData2 + facet_wrap(~cohort)
twinData3 + facet_wrap(~zygosity)## Warning: Removed 141 rows containing missing values (geom_point).
twinData4 <- twinData3 + facet_wrap(~zygosity)
twinData4## Warning: Removed 141 rows containing missing values (geom_point).
## Yes there is also a positive coorelation between the two added
variables cohort and zygosity
library(broom)
twinData %>%
group_by(cohort,zygosity) %>%
do(tidy(cor.test(~ ht1 + ht2, alternative = "greater", data = . ))) ## # A tibble: 10 × 10
## # Groups: cohort, zygosity [10]
## cohort zygosity estimate statistic p.value parame…¹ conf.…² conf.…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr>
## 1 older MZFF 0.859 42.6 1.66e-189 642 0.841 1 Pears…
## 2 older MZMM 0.907 36.4 1.09e-109 286 0.888 1 Pears…
## 3 older DZFF 0.456 10.1 1.06e- 21 387 0.388 1 Pears…
## 4 older DZMM 0.510 6.97 5.84e- 11 138 0.399 1 Pears…
## 5 older DZOS 0.381 8.00 7.65e- 15 378 0.306 1 Pears…
## 6 younger MZFF 0.877 42.7 8.60e-177 547 0.860 1 Pears…
## 7 younger MZMM 0.883 30.0 3.36e- 86 256 0.858 1 Pears…
## 8 younger DZFF 0.440 9.02 7.36e- 18 339 0.365 1 Pears…
## 9 younger DZMM 0.350 5.15 3.32e- 7 190 0.241 1 Pears…
## 10 younger DZOS 0.428 10.4 2.19e- 23 484 0.365 1 Pears…
## # … with 1 more variable: alternative <chr>, and abbreviated variable names
## # ¹​parameter, ²​conf.low, ³​conf.hightable(twinData$zygosity)##
## MZFF MZMM DZFF DZMM DZOS
## 1232 567 751 352 906table(twinData$cohort)##
## older younger
## 1899 1909library(broom)
cor_result <- cor.test(~ ht1 + ht2, alternative = "greater", data = twinData)
tidy_cor_result <- tidy(cor_result)
str(tidy_cor_result)## tibble [1 × 8] (S3: tbl_df/tbl/data.frame)
## $ estimate : Named num 0.619
## ..- attr(*, "names")= chr "cor"
## $ statistic : Named num 47.8
## ..- attr(*, "names")= chr "t"
## $ p.value : num 0
## $ parameter : Named int 3665
## ..- attr(*, "names")= chr "df"
## $ conf.low : num 0.602
## $ conf.high : num 1
## $ method : chr "Pearson's product-moment correlation"
## $ alternative: chr "greater"twinData %>%
group_by(cohort,zygosity) %>%
do(tidy(cor.test(~ ht1 + ht2, alternative = "greater", data = .)))## # A tibble: 10 × 10
## # Groups: cohort, zygosity [10]
## cohort zygosity estimate statistic p.value parame…¹ conf.…² conf.…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr>
## 1 older MZFF 0.859 42.6 1.66e-189 642 0.841 1 Pears…
## 2 older MZMM 0.907 36.4 1.09e-109 286 0.888 1 Pears…
## 3 older DZFF 0.456 10.1 1.06e- 21 387 0.388 1 Pears…
## 4 older DZMM 0.510 6.97 5.84e- 11 138 0.399 1 Pears…
## 5 older DZOS 0.381 8.00 7.65e- 15 378 0.306 1 Pears…
## 6 younger MZFF 0.877 42.7 8.60e-177 547 0.860 1 Pears…
## 7 younger MZMM 0.883 30.0 3.36e- 86 256 0.858 1 Pears…
## 8 younger DZFF 0.440 9.02 7.36e- 18 339 0.365 1 Pears…
## 9 younger DZMM 0.350 5.15 3.32e- 7 190 0.241 1 Pears…
## 10 younger DZOS 0.428 10.4 2.19e- 23 484 0.365 1 Pears…
## # … with 1 more variable: alternative <chr>, and abbreviated variable names
## # ¹​parameter, ²​conf.low, ³​conf.hightwinData %>%
group_by(cohort,zygosity) %>%
do(tidy(cor.test(~ ht1 + ht2, alternative = "greater", data = .)))## # A tibble: 10 × 10
## # Groups: cohort, zygosity [10]
## cohort zygosity estimate statistic p.value parame…¹ conf.…² conf.…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr>
## 1 older MZFF 0.859 42.6 1.66e-189 642 0.841 1 Pears…
## 2 older MZMM 0.907 36.4 1.09e-109 286 0.888 1 Pears…
## 3 older DZFF 0.456 10.1 1.06e- 21 387 0.388 1 Pears…
## 4 older DZMM 0.510 6.97 5.84e- 11 138 0.399 1 Pears…
## 5 older DZOS 0.381 8.00 7.65e- 15 378 0.306 1 Pears…
## 6 younger MZFF 0.877 42.7 8.60e-177 547 0.860 1 Pears…
## 7 younger MZMM 0.883 30.0 3.36e- 86 256 0.858 1 Pears…
## 8 younger DZFF 0.440 9.02 7.36e- 18 339 0.365 1 Pears…
## 9 younger DZMM 0.350 5.15 3.32e- 7 190 0.241 1 Pears…
## 10 younger DZOS 0.428 10.4 2.19e- 23 484 0.365 1 Pears…
## # … with 1 more variable: alternative <chr>, and abbreviated variable names
## # ¹​parameter, ²​conf.low, ³​conf.highsig_twin_cor <- twinData %>%
group_by(cohort,zygosity) %>%
do(tidy( cor.test(~ ht1 + ht2, alternative = "greater" , data = . ))) %>%
arrange(desc(estimate)) %>%
mutate(Greater0.5 = ifelse(estimate>0.5,"Yes","No"))
sig_twin_cor## # A tibble: 10 × 11
## # Groups: cohort, zygosity [10]
## cohort zygosity estimate statistic p.value parame…¹ conf.…² conf.…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr>
## 1 older MZMM 0.907 36.4 1.09e-109 286 0.888 1 Pears…
## 2 younger MZMM 0.883 30.0 3.36e- 86 256 0.858 1 Pears…
## 3 younger MZFF 0.877 42.7 8.60e-177 547 0.860 1 Pears…
## 4 older MZFF 0.859 42.6 1.66e-189 642 0.841 1 Pears…
## 5 older DZMM 0.510 6.97 5.84e- 11 138 0.399 1 Pears…
## 6 older DZFF 0.456 10.1 1.06e- 21 387 0.388 1 Pears…
## 7 younger DZFF 0.440 9.02 7.36e- 18 339 0.365 1 Pears…
## 8 younger DZOS 0.428 10.4 2.19e- 23 484 0.365 1 Pears…
## 9 older DZOS 0.381 8.00 7.65e- 15 378 0.306 1 Pears…
## 10 younger DZMM 0.350 5.15 3.32e- 7 190 0.241 1 Pears…
## # … with 2 more variables: alternative <chr>, Greater0.5 <chr>, and abbreviated
## # variable names ¹​parameter, ²​conf.low, ³​conf.highsort.sig_twin_cor <- sig_twin_cor[order("Yes"),]
sort.sig_twin_cor[1:10,]## # A tibble: 10 × 11
## # Groups: cohort, zygosity [2]
## cohort zygosity estimate statistic p.value parame…¹ conf.…² conf.…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr>
## 1 older MZMM 0.907 36.4 1.09e-109 286 0.888 1 Pears…
## 2 <NA> <NA> NA NA NA NA NA NA <NA>
## 3 <NA> <NA> NA NA NA NA NA NA <NA>
## 4 <NA> <NA> NA NA NA NA NA NA <NA>
## 5 <NA> <NA> NA NA NA NA NA NA <NA>
## 6 <NA> <NA> NA NA NA NA NA NA <NA>
## 7 <NA> <NA> NA NA NA NA NA NA <NA>
## 8 <NA> <NA> NA NA NA NA NA NA <NA>
## 9 <NA> <NA> NA NA NA NA NA NA <NA>
## 10 <NA> <NA> NA NA NA NA NA NA <NA>
## # … with 2 more variables: alternative <chr>, Greater0.5 <chr>, and abbreviated
## # variable names ¹​parameter, ²​conf.low, ³​conf.hightwinData4 <- twinData3 + facet_wrap(~cohort + ~zygosity) + aes(x = ht1, y = ht2, color = cohort)
twinData4 ## Warning: Removed 141 rows containing missing values (geom_point).
hist <- ggplot(twinData, aes(x=ht1, y=ht2)) + geom_boxplot()
hist2 <- hist + facet_wrap(~cohort)
hist3 <- hist2 + facet_wrap(~zygosity)
hist4 <- hist3 + facet_wrap(~cohort + ~zygosity)
hist4## Warning: Continuous x aesthetic -- did you forget aes(group=...)?## Warning: Removed 76 rows containing missing values (stat_boxplot).## Warning: Removed 65 rows containing non-finite values (stat_boxplot).
twin_t <-
twinData %>%
select(cohort,zygosity,ht1,ht2) %>%
group_by(cohort,zygosity) %>%
do(tidy(t.test(.$ht1, .$ht2, paired = TRUE)))
twin_t## # A tibble: 10 × 10
## # Groups: cohort, zygosity [10]
## cohort zygosity estimate stati…¹ p.value param…² conf.low conf.h…³ method
## <chr> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <chr>
## 1 older MZFF 0.00127 0.953 3.41e- 1 643 -1.35e-3 0.00389 Paire…
## 2 older MZMM -0.000239 -0.131 8.95e- 1 287 -3.82e-3 0.00334 Paire…
## 3 older DZFF 0.00314 0.916 3.60e- 1 388 -3.60e-3 0.00988 Paire…
## 4 older DZMM -0.00635 -1.20 2.31e- 1 139 -1.68e-2 0.00409 Paire…
## 5 older DZOS -0.141 -40.5 9.99e-140 379 -1.48e-1 -0.134 Paire…
## 6 younger MZFF 0.000180 0.128 8.98e- 1 548 -2.58e-3 0.00294 Paire…
## 7 younger MZMM 0.00128 0.635 5.26e- 1 257 -2.70e-3 0.00526 Paire…
## 8 younger DZFF 0.00761 1.94 5.30e- 2 340 -9.86e-5 0.0153 Paire…
## 9 younger DZMM 0.00213 0.376 7.07e- 1 191 -9.03e-3 0.0133 Paire…
## 10 younger DZOS -0.143 -43.2 3.16e-168 485 -1.49e-1 -0.136 Paire…
## # … with 1 more variable: alternative <chr>, and abbreviated variable names
## # ¹​statistic, ²​parameter, ³​conf.highhist4 + geom_boxplot(color = "red", fill = "orange", alpha = 0.05)## Warning: Continuous x aesthetic -- did you forget aes(group=...)?## Warning: Removed 76 rows containing missing values (stat_boxplot).## Warning: Removed 65 rows containing non-finite values (stat_boxplot).## Warning: Continuous x aesthetic -- did you forget aes(group=...)?## Warning: Removed 76 rows containing missing values (stat_boxplot).## Warning: Removed 65 rows containing non-finite values (stat_boxplot).
hist4## Warning: Continuous x aesthetic -- did you forget aes(group=...)?## Warning: Removed 76 rows containing missing values (stat_boxplot).## Warning: Removed 65 rows containing non-finite values (stat_boxplot).
hist_9 <- ggplot(twin_t, aes(x=p.value, y=parameter)) + geom_boxplot()
hist_9 + facet_wrap(~cohort + ~zygosity)## Warning: Continuous x aesthetic -- did you forget aes(group=...)?
hist_9## Warning: Continuous x aesthetic -- did you forget aes(group=...)?