Libraries
library(tidyverse)
## Warning in system("timedatectl", intern = TRUE): running command 'timedatectl'
## had status 1
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.6 ✓ dplyr 1.0.7
## ✓ tidyr 1.1.4 ✓ stringr 1.4.0
## ✓ readr 2.1.1 ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(plotly)
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
library(ggbeeswarm)
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
## The following object is masked from 'package:purrr':
##
## some
library(rstatix)
##
## Attaching package: 'rstatix'
## The following object is masked from 'package:stats':
##
## filter
library(SimComp)
library(readxl)
library(dplyr)
library(pwr)
Dataset
dataset_exp_0_bm_imprint <- read_excel (path = "/mnt/c/Users/Toshiya Matsushima/OneDrive/R projects/BM_project/data analysis using Rstudio/BM_project_dataset.xlsx", sheet = "exp_0_preliminary_red_yellow")
dataset_exp_0_bm_imprint
## # A tibble: 42 × 11
## id batch run bm imprint sex drug dose label colour dose_level
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <chr> <chr> <chr> <chr> <dbl>
## 1 id5 1 4076 137 516 1 00_ctrl 0_null 00_re… 0_red 0
## 2 id6 1 2079 54 134 0 00_ctrl 0_null 00_re… 0_red 0
## 3 id7 1 3701 361 578 1 00_ctrl 0_null 00_re… 0_red 0
## 4 id8 1 615 113 589 1 00_ctrl 0_null 00_re… 0_red 0
## 5 id10 4 602 262 584 1 00_ctrl 0_null 00_re… 0_red 0
## 6 id13 4 3073 288 584 0 00_ctrl 0_null 00_re… 0_red 0
## 7 id17 4 511 121 -12 1 00_ctrl 0_null 00_re… 0_red 0
## 8 id19 4 4610 198 437 0 00_ctrl 0_null 00_re… 0_red 0
## 9 id22 5 748 159 591 0 00_ctrl 0_null 00_re… 0_red 0
## 10 id23 5 32 159 447 0 00_ctrl 0_null 00_re… 0_red 0
## # … with 32 more rows
dataset_exp_0_bm_imprint %>%
group_by (label) %>%
get_summary_stats(bm)
## # A tibble: 4 × 14
## label variable n min max median q1 q3 iqr mad mean sd
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 00_red_… bm 10 54 361 159 125 246 121 63.0 185. 93.0
## 2 01_red_… bm 11 -292 564 129 -37 314. 352. 301. 129. 274.
## 3 02_yell… bm 10 -189 567 444. 97.5 532. 435 162. 325. 284.
## 4 03_yell… bm 11 -366 504 167 88 343 255 219. 165. 250.
## # … with 2 more variables: se <dbl>, ci <dbl>
dataset_exp_0_bm_imprint %>%
group_by (label) %>%
get_summary_stats(imprint)
## # A tibble: 4 × 14
## label variable n min max median q1 q3 iqr mad mean sd
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 00_red_… imprint 10 -12 591 547 440. 584 144. 63.8 445. 213.
## 2 01_red_… imprint 11 -311 596 146 -94 378 472 368. 130 292.
## 3 02_yell… imprint 10 -37 522 188 -9 304. 314. 288. 192. 211.
## 4 03_yell… imprint 11 -597 592 -262 -569 0.5 570. 430. -185. 414.
## # … with 2 more variables: se <dbl>, ci <dbl>
x-y plot in all data
fig <- ggplot(data = dataset_exp_0_bm_imprint, mapping = aes(x=bm, y=imprint, color=label))+
geom_point(size=3)+
geom_vline(xintercept=0, linetype = "dotted")+
geom_hline(yintercept=0, linetype = "dotted")+
geom_smooth(method = "lm")+
xlim(-600, 600)+
ylim(-600, 600)+
theme_classic()
fig
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 15 rows containing missing values (geom_smooth).
ggsave(plot = fig, filename = "exp_0_imprint_vs_bm_control_vpa.png", dpi = 300, height = 10, width = 13, units = "cm")
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 15 rows containing missing values (geom_smooth).
imprinting of all data, linear fitting
fit_exp_0_imprint_bm <- lm (imprint ~ bm * colour, data = dataset_exp_0_bm_imprint)
summary (fit_exp_0_imprint_bm)
##
## Call:
## lm(formula = imprint ~ bm * colour, data = dataset_exp_0_bm_imprint)
##
## Residuals:
## Min 1Q Median 3Q Max
## -641.61 -238.01 23.26 211.82 532.18
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 149.3820 89.1486 1.676 0.1020
## bm 0.8369 0.3506 2.387 0.0221 *
## colour1_yellow -215.2363 130.3292 -1.651 0.1069
## bm:colour1_yellow -0.5876 0.4394 -1.337 0.1891
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 322.7 on 38 degrees of freedom
## Multiple R-squared: 0.2805, Adjusted R-squared: 0.2237
## F-statistic: 4.937 on 3 and 38 DF, p-value: 0.005418
car::Anova(fit_exp_0_imprint_bm)
## Anova Table (Type II tests)
##
## Response: imprint
## Sum Sq Df F value Pr(>F)
## bm 499297 1 4.7960 0.034735 *
## colour 1074371 1 10.3198 0.002681 **
## bm:colour 186161 1 1.7882 0.189100
## Residuals 3956096 38
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fig <- ggplot(data = dataset_exp_0_bm_imprint, mapping = aes(x=bm, y=imprint, color=label))+
geom_point(size=3)+
facet_wrap(~drug)+
geom_vline(xintercept=0, linetype = "dotted")+
geom_hline(yintercept=0, linetype = "dotted")+
geom_smooth(method = "lm")+
xlim(-400, 600)+
ylim(-400, 600)+
theme_classic()
fig
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 4 rows containing non-finite values (stat_smooth).
## Warning: Removed 4 rows containing missing values (geom_point).
## Warning: Removed 15 rows containing missing values (geom_smooth).
ggsave(plot = fig, filename = "exp_2_imprint_vs_bm_all_wrap.png", dpi = 300, height = 10, width = 15, units = "cm")
## `geom_smooth()` using formula 'y ~ x'
## Warning: Removed 4 rows containing non-finite values (stat_smooth).
## Warning: Removed 4 rows containing missing values (geom_point).
## Warning: Removed 15 rows containing missing values (geom_smooth).
Analysis of BM
BM all data
fig <- ggplot(data = dataset_exp_0_bm_imprint, mapping = aes(x=label, y=bm))+
geom_hline(yintercept=0, linetype = "dotted")+
geom_boxplot()+
ylim(-400, 600)+
geom_quasirandom(shape=16, size=3, colour="black")+
theme_classic()+
theme(axis.text.x = element_text(angle = 45, hjust = 1))
fig
ggsave(plot = fig, filename = "exp_0_bm_label_all.png", dpi = 300, height = 10, width = 7.5, units = "cm")
Two-way ANOVA of BM
m <- lm(bm ~ drug * colour, data = dataset_exp_0_bm_imprint)
anova(m)
## Analysis of Variance Table
##
## Response: bm
## Df Sum Sq Mean Sq F value Pr(>F)
## drug 1 121363 121363 2.1140 0.1542
## colour 1 76288 76288 1.3289 0.2562
## drug:colour 1 28105 28105 0.4896 0.4884
## Residuals 38 2181492 57408
Analysis of imprinting
imprinting all data
fig <- ggplot(data = dataset_exp_0_bm_imprint, mapping = aes(x=label, y=imprint))+
geom_hline(yintercept=0, linetype = "dotted")+
geom_boxplot()+
ylim(-400, 600)+
geom_quasirandom(shape=16, size=3, colour="black")+
theme_classic()+
theme(axis.text.x = element_text(angle = 45, hjust = 1))
fig
## Warning: Removed 4 rows containing non-finite values (stat_boxplot).
## Warning: Removed 4 rows containing missing values (position_quasirandom).
ggsave(plot = fig, filename = "exp_0_imprint_label_all.png", dpi = 300, height = 10, width = 7.5, units = "cm")
## Warning: Removed 4 rows containing non-finite values (stat_boxplot).
## Warning: Removed 4 rows containing missing values (position_quasirandom).
Two-way ANOVA of imprinting
m <- lm(imprint ~ drug * colour, data = dataset_exp_0_bm_imprint)
anova(m)
## Analysis of Variance Table
##
## Response: imprint
## Df Sum Sq Mean Sq F value Pr(>F)
## drug 1 1252323 1252323 14.0827 0.0005838 ***
## colour 1 856572 856572 9.6324 0.0036004 **
## drug:colour 1 10032 10032 0.1128 0.7388093
## Residuals 38 3379199 88926
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1