DataSet and R Code for Differential Astrocyte Activation by Taenia solium Antigens: Specific Induction of IL-10 and IL-1β by Excretory/Secretory (ES) Products

Statistical analysis Gonzales et. al. Dr. Guillermo Fernandez 2025-12-16 Contents Analysis of cell type distribution (Figure 1) 1 Cytokines gene expressión analysis (Figures 2 and 3) 5 DIV-5analysis ............................................... 6 IL-10.................................................. 6 IL-1B ................................................. 11 TGF-B................................................. 15 DIV-14analysis............................................... 20 IL-10.................................................. 20 IL-1B ................................................. 25 TGF-B................................................. 29 Analysis of cell type distribution (Figure 1) setwd("/home/mi_metodo/Documentos/Laboratorio/UPSJB/Colab/Paper Kevin/estadistica dic 2025/") dat <- read.table("data_cell_type.csv",sep= ",",dec=".",header = TRUE) str(dat) ## 'data.frame': 100 obs. of 5 variables: ## $ DIV : int 14 14 14 14 14 14 14 14 14 14 ... ## $ procedure: chr "a" "a" "a" "a" ... ## $ field : int 1 2 3 4 5 6 7 8 9 10 ... ## $ cell.type: chr "astrocyte" "astrocyte" "astrocyte" "astrocyte" ... ## $ value : int 411 891 671 393 820 425 428 427 445 281 ... library(plyr) library(dplyr) library(Rmisc) library(sciplot) 1 library(ggpubr) library(ggplot2) library(FSA) library (multcompView) library(patchwork) library(gridExtra) library(coin) dat$DIV <- as.factor(dat$DIV) dat$cell.type <- as.factor(dat$cell.type) data1<- dat %>% filter(DIV =='5') data1 ## DIV procedure field cell.type value ## 1 5 a 1 astrocyte 631 ## 2 5 a 2 astrocyte 547 ## 3 5 a 3 astrocyte 681 ## 4 5 a 4 astrocyte 602 ## 5 5 a 5 astrocyte 812 ## 6 5 a 6 astrocyte 809 ## 7 5 a 7 astrocyte 613 ## 8 5 a 8 astrocyte 843 ## 9 5 a 9 astrocyte 823 ## 10 5 a 10 astrocyte 787 ## 11 5 b 1 astrocyte 917 ## 12 5 b 2 astrocyte 946 ## 13 5 b 3 astrocyte 700 ## 14 5 b 4 astrocyte 734 ## 15 5 b 5 astrocyte 802 ## 16 5 b 6 astrocyte 833 ## 17 5 b 7 astrocyte 809 ## 18 5 b 8 astrocyte 662 ## 19 5 b 9 astrocyte 679 ## 20 5 b 10 astrocyte 606 ## 21 5 a 1 neuron 240 ## 22 5 a 2 neuron 201 ## 23 5 a 3 neuron 186 ## 24 5 a 4 neuron 381 ## 25 5 a 5 neuron 284 ## 26 5 a 6 neuron 170 ## 27 5 a 7 neuron 351 ## 28 5 a 8 neuron 265 ## 29 5 a 9 neuron 293 ## 30 5 a 10 neuron 254 ## 31 5 b 1 neuron 208 ## 32 5 b 2 neuron 268 ## 33 5 b 3 neuron 186 ## 34 5 b 4 neuron 227 ## 35 5 b 5 neuron 300 ## 36 5 b 6 neuron 300 ## 37 5 b 7 neuron 387 ## 38 5 b 8 neuron 250 2 ## 39 5 b 9 neuron 274 ## 40 5 b 10 neuron 282 ## 41 5 c 1 astrocyte 785 ## 42 5 c 2 astrocyte 505 ## 43 5 c 3 astrocyte 933 ## 44 5 c 4 astrocyte 878 ## 45 5 c 5 astrocyte 617 ## 46 5 c 6 astrocyte 634 ## 47 5 c 7 astrocyte 837 ## 48 5 c 8 astrocyte 726 ## 49 5 c 9 astrocyte 703 ## 50 5 c 10 astrocyte 762 ## 51 5 c 1 neuron 327 ## 52 5 c 2 neuron 195 ## 53 5 c 3 neuron 435 ## 54 5 c 4 neuron 285 ## 55 5 c 5 neuron 262 ## 56 5 c 6 neuron 171 ## 57 5 c 7 neuron 197 ## 58 5 c 8 neuron 217 ## 59 5 c 9 neuron 248 ## 60 5 c 10 neuron 175 # Creation of the linear model aa <- lm(value~cell.type, data = data1) modelo <- aov(aa) summary(modelo) ## Df Sum Sq Mean Sq F value Pr(>F) ## cell.type 1 3454560 3454560 384 <2e-16 *** ## Residuals 58 521820 8997 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # Normality and homocedasticity test # Diagnostic plots layout(matrix(c(1,2,3,4),2,2)) plot(aa) 3 300 400 500 600 700 −200 100 Fitted values Residuals Residuals vs Fitted 42 12 2 −2 −1 0 1 2 −2 0 2 Theoretical Quantiles Standardized residuals Q−Q Residuals 42 12 2 300 400 500 600 700 0.0 1.0 Fitted values Standardized residuals Scale−Location 42 12 2 −3 0 2 Factor Level Combinations Standardized residuals astrocyte neuron cell.type : Constant Leverage: Residuals vs Factor Levels 42 12 2 # Shapiro Test for normality r=residuals (modelo) shapiro.test(r) ## ## Shapiro-Wilk normality test ## ## data: r ## W = 0.9923, p-value = 0.9698 # Homocedasticity test bartlett.test(value~cell.type,data=data1) ## ## Bartlett test of homogeneity of variances ## ## data: value by cell.type ## Bartlett's K-squared = 7.8916, df = 1, p-value = 0.004967 # Mann-Whitney test wilcox_test(value ~cell.type, data = data1, distribution = "exact",conf.int=0.95) 4 ## ## Exact Wilcoxon-Mann-Whitney Test ## ## data: value by cell.type (astrocyte, neuron) ## Z = 6.6533, p-value < 2.2e-16 ## alternative hypothesis: true mu is not equal to 0 ## 95 percent confidence interval: ## 427 539 ## sample estimates: ## difference in location ## 485 # Descriptive statistics resumen <- data1 %>% dplyr::filter(cell.type %in% c("astrocyte","neuron")) %>% dplyr::group_by(cell.type) %>% dplyr::summarise( n = sum(!is.na(value)), mediana = median(value, na.rm = TRUE), IQR = IQR(value, na.rm = TRUE), p25 = quantile(value, 0.25,na.rm = TRUE), p75 = quantile(value, 0.75,na.rm = TRUE), min = min(value, na.rm = TRUE), max = max(value, na.rm = TRUE), .groups = "drop" ) print(resumen) ## # A tibble: 2 x 8 ## cell.type n mediana IQR p25 p75 min max ## <fct> <int> <dbl> <dbl> <dbl> <dbl> <int> <int> ## 1 astrocyte 30 748 179. 641 820. 505 946 ## 2 neuron 30 258 88.2 203. 291 170 435 Cytokines gene expressión analysis (Figures 2 and 3) dat <- read.table("data_cytokine.csv",sep= ",",dec=".",header = TRUE) dat$div <- as.factor(dat$div) dat$treatment <- as.factor(dat$treatment) dat$n<- as.character(dat$n) dat$cytokine <- as.factor(dat$cytokine) str(dat) ## 'data.frame': 216 obs. of 5 variables: ## $ div : Factor w/ 2 levels "5","14": 2 2 2 2 2 2 2 2 2 2 ... ## $ treatment: Factor w/ 4 levels "ctrl","es","lps",..: 1 1 1 1 1 1 1 1 1 1 ... 5 ## $ n : chr "1" "2" "3" "4" ... ## $ cytokine : Factor w/ 3 levels "il-10","il-1b",..: 1 1 1 1 1 1 1 1 1 2 ... ## $ value : num 11.6 10.1 11.2 14.2 13.9 ... sapply(dat, function(x) if(is.factor(x)) levels(x)) ## $div ## [1] "5" "14" ## ## $treatment ## [1] "ctrl" "es" "lps" "total" ## ## $n ## NULL ## ## $cytokine ## [1] "il-10" "il-1b" "tgfb" ## ## $value ## NULL DIV-5 analysis IL-10 # Filter DIV5 and IL-10 data1<- dat %>% filter(div =='5')%>% filter(cytokine =='il-10') # Obtained dataframe after filtering data1 ## div treatment n cytokine value ## 1 5 ctrl 1 il-10 17.30 ## 2 5 ctrl 2 il-10 16.68 ## 3 5 ctrl 3 il-10 17.69 ## 4 5 ctrl 4 il-10 14.13 ## 5 5 ctrl 5 il-10 11.18 ## 6 5 ctrl 6 il-10 12.62 ## 7 5 ctrl 7 il-10 13.32 ## 8 5 ctrl 8 il-10 13.67 ## 9 5 ctrl 9 il-10 14.37 ## 10 5 lps 1 il-10 10.67 ## 11 5 lps 2 il-10 10.49 ## 12 5 lps 3 il-10 10.07 ## 13 5 lps 4 il-10 8.63 ## 14 5 lps 5 il-10 8.03 ## 15 5 lps 6 il-10 8.65 ## 16 5 lps 7 il-10 8.70 ## 17 5 lps 8 il-10 9.01 6 ## 18 5 lps 9 il-10 8.51 ## 19 5 es 1 il-10 10.49 ## 20 5 es 2 il-10 10.86 ## 21 5 es 3 il-10 10.69 ## 22 5 es 4 il-10 11.43 ## 23 5 es 5 il-10 9.28 ## 24 5 es 6 il-10 12.06 ## 25 5 es 7 il-10 9.47 ## 26 5 es 8 il-10 9.21 ## 27 5 es 9 il-10 9.19 ## 28 5 total 1 il-10 16.26 ## 29 5 total 2 il-10 17.00 ## 30 5 total 3 il-10 17.58 ## 31 5 total 4 il-10 12.60 ## 32 5 total 5 il-10 12.49 ## 33 5 total 6 il-10 12.87 ## 34 5 total 7 il-10 10.19 ## 35 5 total 8 il-10 10.20 ## 36 5 total 9 il-10 11.29 # Creation of the linear model aa <- lm(value~treatment, data = data1) modelo <- aov(aa) summary(modelo) ## Df Sum Sq Mean Sq F value Pr(>F) ## treatment 3 172.0 57.34 15.13 2.63e-06 *** ## Residuals 32 121.3 3.79 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # Normality and homocedasticity test # Diagnostic plots layout(matrix(c(1,2,3,4),2,2)) # optional layout plot(aa) # diagnostic plots 7 9 10 11 12 13 14 −4 0 4 Fitted values Residuals Residuals vs Fitted 30 29 5 −2 −1 0 1 2 −2 0 2 Theoretical Quantiles Standardized residuals Q−Q Residuals 30 29 5 9 10 11 12 13 14 0.0 1.0 Fitted values Standardized residuals Scale−Location 30 29 5 −2 0 2 Factor Level Combinations Standardized residuals ctrl es lps total treatment : Constant Leverage: Residuals vs Factor Levels 30 29 5 # Shapiro Test for normality r=residuals (modelo) shapiro.test(r) ## ## Shapiro-Wilk normality test ## ## data: r ## W = 0.9551, p-value = 0.1513 # Homocedasticity test bartlett.test(value~treatment,data=data1) ## ## Bartlett test of homogeneity of variances ## ## data: value by treatment ## Bartlett's K-squared = 12.136, df = 3, p-value = 0.006932 # Descriptive statistics 8 # Mean plus SE bb<- summarySE(data1, measurevar = "value",groupvars = "treatment") bb ## treatment N value sd se ci ## 1 ctrl 9 14.551111 2.2230859 0.7410286 1.7088151 ## 2 es 9 10.297778 1.0621767 0.3540589 0.8164613 ## 3 lps 9 9.195556 0.9578767 0.3192922 0.7362892 ## 4 total 9 13.386667 2.8582250 0.9527417 2.1970262 # median and IQR ranges resumen_dplyr <- data1 %>% group_by(treatment) %>% summarise( n = sum(!is.na(value)), mediana = median(value, na.rm = TRUE), p25 = quantile(value, 0.25,na.rm = TRUE), p75 = quantile(value, 0.75,na.rm = TRUE), IQR = IQR(value, na.rm = TRUE), .groups = "drop" ) resumen_dplyr <- resumen_dplyr %>% mutate(mediana_IQR = sprintf("%.2f [%.2f–%.2f]", mediana, p25, p75)) resumen_dplyr ## # A tibble: 4 x 7 ## treatment n mediana p25 p75 IQR mediana_IQR ## <fct> <int> <dbl> <dbl> <dbl> <dbl> <chr> ## 1 ctrl 9 14.1 13.3 16.7 3.36 14.13 [13.32–16.68] ## 2 es 9 10.5 9.28 10.9 1.58 10.49 [9.28–10.86] ## 3 lps 9 8.7 8.63 10.1 1.44 8.70 [8.63–10.07] ## 4 total 9 12.6 11.3 16.3 4.97 12.60 [11.29–16.26] # Boxplot ggplot(data1, aes(x = treatment, y = value, fill = treatment)) + geom_boxplot(outlier.alpha = 0.7)+ scale_fill_manual(values = c(ctrl = "#05090AFF",es = "#833437FF",lps = "#67B9E9FF", total = "#E9D097FF")) + labs(title = "DIV 5Il-10", x = "Treatment",y = "Delta CT")+ theme_minimal() 9 ## 11 5 lps 2 tgfb 8.32 ## 12 5 lps 3 tgfb 7.94 ## 13 5 lps 4 tgfb 6.48 ## 14 5 lps 5 tgfb 7.07 ## 15 5 lps 6 tgfb 6.55 ## 16 5 lps 7 tgfb 6.87 ## 17 5 lps 8 tgfb 7.35 ## 18 5 lps 9 tgfb 6.94 ## 19 5 es 1 tgfb 7.69 ## 20 5 es 2 tgfb 7.72 ## 21 5 es 3 tgfb 7.67 ## 22 5 es 4 tgfb 7.56 ## 23 5 es 5 tgfb 7.53 ## 24 5 es 6 tgfb 9.63 ## 25 5 es 7 tgfb 7.36 ## 26 5 es 8 tgfb 6.87 ## 27 5 es 9 tgfb 7.02 ## 28 5 total 1 tgfb 7.10 ## 29 5 total 2 tgfb 6.93 ## 30 5 total 3 tgfb 7.25 ## 31 5 total 4 tgfb 7.19 ## 32 5 total 5 tgfb 6.27 ## 33 5 total 6 tgfb 7.65 ## 34 5 total 7 tgfb 5.76 ## 35 5 total 8 tgfb 5.70 ## 36 5 total 9 tgfb 6.10 # Create the linear model aa <- lm(value~treatment, data = data1) modelo <- aov(aa) summary(modelo) ## Df Sum Sq Mean Sq F value Pr(>F) ## treatment 3 5.76 1.920 4.033 0.0153 * ## Residuals 32 15.23 0.476 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # Normality and homocedasticity test # Diagnostic plots layout(matrix(c(1,2,3,4),2,2)) # optional layout plot(aa) # diagnostic plots 16 6.8 7.0 7.2 7.4 7.6 −1.0 1.0 Fitted values Residuals Residuals vs Fitted 24 5 11 −2 −1 0 1 2 −1 1 3 Theoretical Quantiles Standardized residuals Q−Q Residuals 24 5 11 6.8 7.0 7.2 7.4 7.6 0.0 1.0 Fitted values Standardized residuals Scale−Location 24 511 −2 0 2 Factor Level Combinations Standardized residuals ctrl es lps total treatment : Constant Leverage: Residuals vs Factor Levels 24 5 11 # Shapiro Test for normality r=residuals (modelo) shapiro.test(r) ## ## Shapiro-Wilk normality test ## ## data: r ## W = 0.96212, p-value = 0.2497 # Homocedasticity test bartlett.test(value~treatment,data=data1) ## ## Bartlett test of homogeneity of variances ## ## data: value by treatment ## Bartlett's K-squared = 0.99335, df = 3, p-value = 0.8029 # Descriptive statistics 17 # Mean with SE bb<- summarySE(data1, measurevar = "value",groupvars = "treatment") bb ## treatment N value sd se ci ## 1 ctrl 9 6.810000 0.5530597 0.1843532 0.4251193 ## 2 es 9 7.672222 0.7933123 0.2644374 0.6097938 ## 3 lps 9 7.290000 0.6780487 0.2260162 0.5211943 ## 4 total 9 6.661111 0.7135202 0.2378401 0.5484602 # median and IQR ranges resumen_dplyr <- data1 %>% group_by(treatment) %>% summarise( n = sum(!is.na(value)), mediana = median(value, na.rm = TRUE), p25 = quantile(value, 0.25,na.rm = TRUE), p75 = quantile(value, 0.75,na.rm = TRUE), IQR = IQR(value, na.rm = TRUE), .groups = "drop" ) resumen_dplyr <- resumen_dplyr %>% mutate(mediana_IQR = sprintf("%.2f [%.2f–%.2f]", mediana, p25, p75)) resumen_dplyr ## # A tibble: 4 x 7 ## treatment n mediana p25 p75 IQR mediana_IQR ## <fct> <int> <dbl> <dbl> <dbl> <dbl> <chr> ## 1 ctrl 9 6.68 6.62 7.34 0.72 6.68 [6.62–7.34] ## 2 es 9 7.56 7.36 7.69 0.33 7.56 [7.36–7.69] ## 3 lps 9 7.07 6.87 7.94 1.07 7.07 [6.87–7.94] ## 4 total 9 6.93 6.1 7.19 1.09 6.93 [6.10–7.19] # Boxplot ggplot(data1, aes(x = treatment, y = value, fill = treatment)) + geom_boxplot(outlier.alpha = 0.7)+ scale_fill_manual(values = c(ctrl = "#05090AFF",es = "#833437FF",lps = "#67B9E9FF", total = "#E9D097FF")) + labs(title = "DIV 5TGF-B", x = "Treatment",y = "Delta CT")+ theme_minimal() 18 6 7 8 9 ctrl es lps total Treatment Delta CT treatment ctrl es lps total DIV 5− TGF−B # Non -parametric test kw <- kruskal.test(value ~treatment, data = data1) kw ## ## Kruskal-Wallis rank sum test ## ## data: value by treatment ## Kruskal-Wallis chi-squared = 9.304, df = 3, p-value = 0.02551 # post-hoc comparisonDunn test with holm correction dunn <- dunnTest(value ~treatment, data = data1, method = "holm") # Dunn test with 4 decs and asterisc in significant results results <- transform( dunn$res, P.unadj = sprintf("%.4f", P.unadj), P.adj = sprintf("%.4f", P.adj) ) results$signif <- with(dunn$res, ifelse(P.adj <0.001,"***", ifelse(P.adj <0.01,"**", ifelse(P.adj <0.05,"*","ns")))) 19 print(results) ## Comparison Z P.unadj P.adj signif ## 1 ctrl - es -2.4275039 0.0152 0.0760 ns ## 2 ctrl - lps -1.2864652 0.1983 0.5948 ns ## 3 es - lps 1.1410387 0.2539 0.5077 ns ## 4 ctrl - total 0.2684797 0.7883 0.7883 ns ## 5 es - total 2.6959836 0.0070 0.0421 * ## 6 lps - total 1.5549449 0.1200 0.4798 ns DIV-14 analysis IL-10 # Filter DIV14 and IL-10 data1<- dat %>% filter(div =='14')%>% filter(cytokine =='il-10') # Obtained dataframe after filtering data1 ## div treatment n cytokine value ## 1 14 ctrl 1 il-10 11.63 ## 2 14 ctrl 2 il-10 10.08 ## 3 14 ctrl 3 il-10 11.15 ## 4 14 ctrl 4 il-10 14.23 ## 5 14 ctrl 5 il-10 13.88 ## 6 14 ctrl 6 il-10 12.69 ## 7 14 ctrl 7 il-10 11.05 ## 8 14 ctrl 8 il-10 17.03 ## 9 14 ctrl 9 il-10 14.46 ## 10 14 lps 1 il-10 10.17 ## 11 14 lps 2 il-10 9.96 ## 12 14 lps 3 il-10 9.47 ## 13 14 lps 4 il-10 9.13 ## 14 14 lps 5 il-10 9.18 ## 15 14 lps 6 il-10 8.86 ## 16 14 lps 7 il-10 7.21 ## 17 14 lps 8 il-10 8.59 ## 18 14 lps 9 il-10 7.90 ## 19 14 es 1 il-10 8.00 ## 20 14 es 2 il-10 9.73 ## 21 14 es 3 il-10 10.38 ## 22 14 es 4 il-10 10.18 ## 23 14 es 5 il-10 10.86 ## 24 14 es 6 il-10 9.80 ## 25 14 es 7 il-10 6.95 ## 26 14 es 8 il-10 7.54 ## 27 14 es 9 il-10 6.64 20 ## 28 14 total 1 il-10 7.03 ## 29 14 total 2 il-10 9.16 ## 30 14 total 3 il-10 8.36 ## 31 14 total 4 il-10 10.21 ## 32 14 total 5 il-10 10.26 ## 33 14 total 6 il-10 10.98 ## 34 14 total 7 il-10 11.02 ## 35 14 total 8 il-10 11.10 ## 36 14 total 9 il-10 11.10 # Creation of the linear model aa <- lm(value~treatment, data = data1) modelo <- aov(aa) summary(modelo) ## Df Sum Sq Mean Sq F value Pr(>F) ## treatment 3 96.37 32.12 12.41 1.51e-05 *** ## Residuals 32 82.86 2.59 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # Normality and homocedasticity test # Diagnostic plots layout(matrix(c(1,2,3,4),2,2)) # optional layout plot(aa) # diagnostic plots 21 9 10 11 12 13 −2 2 Fitted values Residuals Residuals vs Fitted 8 28 2 −2 −1 0 1 2 −2 0 2 Theoretical Quantiles Standardized residuals Q−Q Residuals 8 28 2 9 10 11 12 13 0.0 1.0 Fitted values Standardized residuals Scale−Location 8 28 2 −2 0 2 Factor Level Combinations Standardized residuals ctrl es lps total treatment : Constant Leverage: Residuals vs Factor Levels 8 28 2 # Shapiro Test for normality r=residuals (modelo) shapiro.test(r) ## ## Shapiro-Wilk normality test ## ## data: r ## W = 0.95306, p-value = 0.1306 # Homocedasticity test bartlett.test(value~treatment,data=data1) ## ## Bartlett test of homogeneity of variances ## ## data: value by treatment ## Bartlett's K-squared = 5.0724, df = 3, p-value = 0.1666 # Descriptive statistics # Mean plus SE 22 bb<- summarySE(data1, measurevar = "value",groupvars = "treatment") bb ## treatment N value sd se ci ## 1 ctrl 9 12.911111 2.1884095 0.7294698 1.6821604 ## 2 es 9 8.897778 1.6100837 0.5366946 1.2376199 ## 3 lps 9 8.941111 0.9441457 0.3147152 0.7257346 ## 4 total 9 9.913333 1.4439962 0.4813321 1.1099537 # median and IQR ranges resumen_dplyr <- data1 %>% group_by(treatment) %>% summarise( n = sum(!is.na(value)), mediana = median(value, na.rm = TRUE), p25 = quantile(value, 0.25,na.rm = TRUE), p75 = quantile(value, 0.75,na.rm = TRUE), IQR = IQR(value, na.rm = TRUE), .groups = "drop" ) resumen_dplyr <- resumen_dplyr %>% mutate(mediana_IQR = sprintf("%.2f [%.2f–%.2f]", mediana, p25, p75)) resumen_dplyr ## # A tibble: 4 x 7 ## treatment n mediana p25 p75 IQR mediana_IQR ## <fct> <int> <dbl> <dbl> <dbl> <dbl> <chr> ## 1 ctrl 9 12.7 11.2 14.2 3.08 12.69 [11.15–14.23] ## 2 es 9 9.73 7.54 10.2 2.64 9.73 [7.54–10.18] ## 3 lps 9 9.13 8.59 9.47 0.880 9.13 [8.59–9.47] ## 4 total 9 10.3 9.16 11.0 1.86 10.26 [9.16–11.02] # Boxplot ggplot(data1, aes(x = treatment, y = value, fill = treatment)) + geom_boxplot(outlier.alpha = 0.7)+ scale_fill_manual(values = c(ctrl = "#05090AFF",es = "#833437FF",lps = "#67B9E9FF", total = "#E9D097FF")) + labs(title = "DIV 14Il-10", x = "Treatment",y = "Delta CT")+ theme_minimal() 23 9 12 15 ctrl es lps total Treatment Delta CT treatment ctrl es lps total DIV 14− Il−10 # Non -parametric test kw <- kruskal.test(value ~treatment, data = data1) kw ## ## Kruskal-Wallis rank sum test ## ## data: value by treatment ## Kruskal-Wallis chi-squared = 18.962, df = 3, p-value = 0.0002783 # post-hoc comparisonDunn test with holm correction dunn <- dunnTest(value ~treatment, data = data1, method = "holm") # Dunn test with 4 decs and asterisc in significant results results <- transform( dunn$res, P.unadj = sprintf("%.4f", P.unadj), P.adj = sprintf("%.4f", P.adj) ) results$signif <- with(dunn$res, ifelse(P.adj <0.001,"***", ifelse(P.adj <0.01,"**", ifelse(P.adj <0.05,"*","ns")))) print(results) 24 ## Comparison Z P.unadj P.adj signif ## 1 ctrl - es 3.6468492 0.0003 0.0013 ** ## 2 ctrl - lps 3.8705823 0.0001 0.0007 *** ## 3 es - lps 0.2237331 0.8230 0.8230 ns ## 4 ctrl - total 2.2820774 0.0225 0.0899 ns ## 5 es - total -1.3647718 0.1723 0.3446 ns ## 6 lps - total -1.5885049 0.1122 0.3365 ns IL-1B # Filter DIV14 and IL-1B data1<- dat %>% filter(div =='14')%>% filter(cytokine =='il-1b') # Obtained dataframe after filtering data1 ## div treatment n cytokine value ## 1 14 ctrl 1 il-1b 4.10 ## 2 14 ctrl 2 il-1b 4.45 ## 3 14 ctrl 3 il-1b 3.52 ## 4 14 ctrl 4 il-1b 9.03 ## 5 14 ctrl 5 il-1b 8.29 ## 6 14 ctrl 6 il-1b 7.39 ## 7 14 ctrl 7 il-1b 7.90 ## 8 14 ctrl 8 il-1b 12.70 ## 9 14 ctrl 9 il-1b 13.67 ## 10 14 lps 1 il-1b 3.49 ## 11 14 lps 2 il-1b 2.26 ## 12 14 lps 3 il-1b 1.63 ## 13 14 lps 4 il-1b 4.94 ## 14 14 lps 5 il-1b 4.42 ## 15 14 lps 6 il-1b 4.63 ## 16 14 lps 7 il-1b 3.93 ## 17 14 lps 8 il-1b 2.73 ## 18 14 lps 9 il-1b 3.49 ## 19 14 es 1 il-1b 3.07 ## 20 14 es 2 il-1b 3.23 ## 21 14 es 3 il-1b 5.91 ## 22 14 es 4 il-1b 5.57 ## 23 14 es 5 il-1b 5.95 ## 24 14 es 6 il-1b 4.82 ## 25 14 es 7 il-1b 2.95 ## 26 14 es 8 il-1b 3.61 ## 27 14 es 9 il-1b 2.69 ## 28 14 total 1 il-1b 1.98 ## 29 14 total 2 il-1b 3.65 ## 30 14 total 3 il-1b 4.04 ## 31 14 total 4 il-1b 5.45 ## 32 14 total 5 il-1b 5.15 25 # Mean with SE bb<- summarySE(data1, measurevar = "value",groupvars = "treatment") bb ## treatment N value sd se ci ## 1 ctrl 9 7.034444 0.8320323 0.2773441 0.6395567 ## 2 es 9 6.763333 0.6088924 0.2029641 0.4680362 ## 3 lps 9 7.186667 0.6225753 0.2075251 0.4785537 ## 4 total 9 6.290000 0.5971181 0.1990394 0.4589856 # median and IQR ranges resumen_dplyr <- data1 %>% group_by(treatment) %>% summarise( n = sum(!is.na(value)), mediana = median(value, na.rm = TRUE), p25 = quantile(value, 0.25,na.rm = TRUE), p75 = quantile(value, 0.75,na.rm = TRUE), IQR = IQR(value, na.rm = TRUE), .groups = "drop" ) resumen_dplyr <- resumen_dplyr %>% mutate(mediana_IQR = sprintf("%.2f [%.2f–%.2f]", mediana, p25, p75)) resumen_dplyr ## # A tibble: 4 x 7 ## treatment n mediana p25 p75 IQR mediana_IQR ## <fct> <int> <dbl> <dbl> <dbl> <dbl> <chr> ## 1 ctrl 9 7.2 6.44 7.4 0.96 7.20 [6.44–7.40] ## 2 es 9 6.88 6.55 7.24 0.690 6.88 [6.55–7.24] ## 3 lps 9 7.07 6.73 7.42 0.690 7.07 [6.73–7.42] ## 4 total 9 6.24 5.67 6.85 1.18 6.24 [5.67–6.85] # Boxplot ggplot(data1, aes(x = treatment, y = value, fill = treatment)) + geom_boxplot(outlier.alpha = 0.7)+ scale_fill_manual(values = c(ctrl = "#05090AFF",es = "#833437FF",lps = "#67B9E9FF", total = "#E9D097FF")) + labs(title = "DIV 14TGF-B", x = "Treatment",y = "Delta CT")+ theme_minimal() 32 6 7 8 ctrl es lps total Treatment Delta CT treatment ctrl es lps total DIV 14− TGF−B # Non -parametric test kw <- kruskal.test(value ~treatment, data = data1) kw ## ## Kruskal-Wallis rank sum test ## ## data: value by treatment ## Kruskal-Wallis chi-squared = 7.1912, df = 3, p-value = 0.06605 # post-hoc comparisonDunn test with holm correction dunn <- dunnTest(value ~treatment, data = data1, method = "holm") # Dunn test with 4 decs and asterisc in significant results results <- transform( dunn$res, P.unadj = sprintf("%.4f", P.unadj), P.adj = sprintf("%.4f", P.adj) ) results$signif <- with(dunn$res, ifelse(P.adj <0.001,"***", ifelse(P.adj <0.01,"**", ifelse(P.adj <0.05,"*","ns")))) print(results) 33 ## Comparison Z P.unadj P.adj signif ## 1 ctrl - es 0.5818933 0.5606 1.0000 ns ## 2 ctrl - lps -0.2573759 0.7969 0.7969 ns ## 3 es - lps -0.8392692 0.4013 1.0000 ns ## 4 ctrl - total 2.1820999 0.0291 0.1455 ns ## 5 es - total 1.6002066 0.1096 0.4382 ns ## 6 lps - total 2.4394757 0.0147 0.0883 ns 34