Following libraries will be useful. If the libraries are not installed, please install them using install.packages()
In a cell culture lab various cellular assays are performed. The package “bioassays” will help to analyse the results of these experiments performed in multiwell plates. The functions in this package can be used to summarise data from any multiwell plate, and by incorporating them in a loop several plates can be analyzed automatically. Two examples are shown in this article. All the csv files used for the examples are availble in the folder exdata (inst/exdata).
To set up a folder as working directory
path1<-tk_choose.dir(getwd(), "Choose the folder for Analysis")
# A window will popup and ask you to select the folder containg the data
setwd(path1)
Read files
Files can be read using the codes below. Files need to be in .csv format. In this examples file names reflect the type of assay used. For example “L_DIFF_P3_72HRS.csv” :L is assay code (Lactate assay), DIFF differentiation (type of cells used), p3: Plate 3 (each plate represent specific compounds tested), 72HRS 72hrs treatment with compound. These information will help to summarise the results.
filelist <-list("L_HEPG2_P3_72HRS.csv","L_HEPG2_P3_24HRS.csv") # list of files
fno<-1 # file number (in fileslist) that is going to be analyzed
result <- data.frame(stringsAsFactors= FALSE) ## An empty dataframe to dump result
zzz <- data.frame(stringsAsFactors= FALSE) ## An empty dataframe to dump result
To read 1st file
filename<-extract_filename(filelist[fno])[1]
filename
#> [1] "L_HEPG2_P3_72HRS.csv"
nickname<-extract_filename(filelist[fno], split="_",end=".csv",remove="",sep="")[2]
nickname
#> [1] "LHEPG2P372HRS"
rawdata<-read.csv(filename,stringsAsFactors = FALSE, strip.white = TRUE,
na.strings = c("NA",""),header = TRUE,skip=1)
head(rawdata)
#> X X1 X2 X3 X4 X5 X6 X7 X8 X9 X10 X11 X12
#> 1 A 0.659 0.649 0.598 0.601 0.541 0.553 0.568 0.519 0.576 0.575 0.583 0.504
#> 2 B 0.442 0.455 0.586 0.563 0.525 0.548 0.511 0.503 0.533 0.559 0.529 0.535
#> 3 C 0.278 0.266 0.491 0.562 0.510 0.473 0.467 0.433 0.382 0.457 0.475 0.510
#> 4 D 0.197 0.199 0.452 0.456 0.421 0.431 0.409 0.401 0.458 0.412 0.408 0.403
#> 5 E 0.177 0.174 0.447 0.437 0.392 0.412 0.368 0.396 0.397 0.358 0.360 0.393
#> 6 F 0.141 0.137 0.277 0.337 0.294 0.279 0.257 0.263 0.262 0.292 0.280 0.300
Reading metadata file
metadata<-read.csv("metafile.csv",stringsAsFactors = FALSE,strip.white = TRUE,
na.strings = c("NA",""),header = TRUE)
head(metadata)
#> row col position type id concentration dilution
#> 1 A 1 A01 STD1 STD 25 NA
#> 2 A 2 A02 STD1 STD 25 NA
#> 3 A 3 A03 S1 Sample NA NA
#> 4 A 4 A04 S1 Sample NA NA
#> 5 A 5 A05 S1 Sample NA NA
#> 6 A 6 A06 S1 Sample NA NA
96 well plates were used for the assay, so it is assumed that in data the rows being from A to G and columns 1 to 12 of a 96 well plate. data2plateformat function is used to label them correctly.
rawdata<-data2plateformat(rawdata,platetype = 96)
head(rawdata)
#> 1 2 3 4 5 6 7 8 9 10 11 12
#> A 0.659 0.649 0.598 0.601 0.541 0.553 0.568 0.519 0.576 0.575 0.583 0.504
#> B 0.442 0.455 0.586 0.563 0.525 0.548 0.511 0.503 0.533 0.559 0.529 0.535
#> C 0.278 0.266 0.491 0.562 0.510 0.473 0.467 0.433 0.382 0.457 0.475 0.510
#> D 0.197 0.199 0.452 0.456 0.421 0.431 0.409 0.401 0.458 0.412 0.408 0.403
#> E 0.177 0.174 0.447 0.437 0.392 0.412 0.368 0.396 0.397 0.358 0.360 0.393
#> F 0.141 0.137 0.277 0.337 0.294 0.279 0.257 0.263 0.262 0.292 0.280 0.300
data2plateformat function uses the column name and row name of the above rawdata to format it as a data frame. So the data should be well formatted before using this function.
An example function is given below which will determine the compound,concentration, type and dilution from the file name
plate_info<-function(file,i){
file<-file[1]
plate<- extract_filename(file,split = "_",end = ".csv", remove = " ", sep=" ")[5]
if(plate == "P2"){
compound<-"CyclosporinA" # Concentration of cyclosporinA used for experiment
concentration<-c(0,1,5,10,15,20,25,50) # Concentration of cyclosporinA used for experiment
type<-c("S1","S2","S3","S4","S5","S6","S7","S8") # sample names of corresponding concentration
dilution<-5
plate_meta<-list(compound=compound,concentration=round(concentration,2),type=type,
dilution=dilution)
}
if(plate == "P3"){
compound<-"Taxol"
concentration<-c(0,0.0125,.025,.05,0.1,1,5,10)
type<-c("S1","S2","S3","S4","S5","S6","S7","S8")
dilution<-5
plate_meta<-list(compound=compound,concentration=round(concentration,2),type=type,
dilution=dilution)
}
if(plate =="p4"){
compound<-c("Cisplatin")
concentration<-c(0,0.5,2,4,8,16,64,"")
type<-c("S1","S2","S3","S4","S5","S6","S7","")
dilution <- 5
plate_meta<-list(compound=compound,concentration=round(concentration,2),type=type,
dilution=dilution)
}
return(plate_meta)
}
plate_details<-plate_info(filelist,1)
plate_details
#> $compound
#> [1] "Taxol"
#>
#> $concentration
#> [1] 0.00 0.01 0.03 0.05 0.10 1.00 5.00 10.00
#>
#> $type
#> [1] "S1" "S2" "S3" "S4" "S5" "S6" "S7" "S8"
#>
#> $dilution
#> [1] 5
These ‘plate details’ can be used to fill metadata using plate_metadata function
metadata1<-plate_metadata(plate_details,metadata,mergeby="type")
head(metadata1)
#> row col type position id dilution concentration compound
#> 1 A 1 STD1 A01 STD NA 25 <NA>
#> 2 A 2 STD1 A02 STD NA 25 <NA>
#> 3 A 3 S1 A03 Sample 5 0 Taxol
#> 4 A 4 S1 A04 Sample 5 0 Taxol
#> 5 A 5 S1 A05 Sample 5 0 Taxol
#> 6 A 6 S1 A06 Sample 5 0 Taxol
For joining metadata and platelayout ‘inner_join’ function is used
data_DF<- dplyr::inner_join(OD_df,metadata1,by=c("row","col","position"))
assign(paste("data",sep="_",nickname),data_DF) # create a copy of data with platename
head(data_DF)
#> row col position value type id dilution concentration compound
#> 1 A 1 A01 0.659 STD1 STD NA 25 <NA>
#> 2 A 2 A02 0.649 STD1 STD NA 25 <NA>
#> 3 A 3 A03 0.598 S1 Sample 5 0 Taxol
#> 4 A 4 A04 0.601 S1 Sample 5 0 Taxol
#> 5 A 5 A05 0.541 S1 Sample 5 0 Taxol
#> 6 A 6 A06 0.553 S1 Sample 5 0 Taxol
Blank can be reduced using ‘reduceblank’ function
data_DF<-reduceblank(data_DF, x_vector =c("All"),blank_vector = c("Blank"), "value")
head(data_DF)
#> row col position value type id dilution concentration compound blankminus
#> 1 A 1 A01 0.659 STD1 STD NA 25 <NA> 0.5545
#> 2 A 2 A02 0.649 STD1 STD NA 25 <NA> 0.5445
#> 3 A 3 A03 0.598 S1 Sample 5 0 Taxol 0.4935
#> 4 A 4 A04 0.601 S1 Sample 5 0 Taxol 0.4965
#> 5 A 5 A05 0.541 S1 Sample 5 0 Taxol 0.4365
#> 6 A 6 A06 0.553 S1 Sample 5 0 Taxol 0.4485
assign(paste("Blkmin",sep="_",nickname),data_DF) # create a copy of data with platename
To filter standards
std<- dplyr::filter(data_DF, data_DF$id=="STD")
std<- aggregate(std$blankminus ~ std$concentration, FUN = mean )
colnames (std) <-c("con", "OD")
head(std)
#> con OD
#> 1 0.39 0.0155
#> 2 0.78 0.0345
#> 3 1.56 0.0710
#> 4 3.13 0.0935
#> 5 6.25 0.1675
#> 6 12.50 0.3440
Calculations for standard curve : nonparametric logistic regression curve
fit1<-nplr::nplr(std$con,std$OD,npars=3,useLog = FALSE)
#npars = 3 for 3 parametric regression curve
#for graph
x1 <- nplr::getX(fit1); y1 <- nplr::getY(fit1)
x2 <- nplr::getXcurve(fit1); y2 <- nplr::getYcurve(fit1)
plot(x1, y1, pch=15, cex=1, col="red", xlab="Concentration",
ylab="Mean OD", main=paste("Standard Curve: ", nickname), cex.main=1)
lines(x2, y2, lwd=3, col="seagreen4")
To evaluate nonparametric logistic regression fitting
params<-nplr::getPar(fit1)$params
nplr::getGoodness(fit1)
#> $gof
#> [1] 0.9919475
#>
#> $wgof
#> [1] 0.9995859
To estimate the values based on logistic regression fitting
estimated_nplr<-estimate(data_DF,colname="blankminus",fitformula=fit1,method="nplr")
#> These values have been replaced by the minimal possible value the model can estimate.
head(estimated_nplr)
#> row col position value type id dilution concentration compound blankminus
#> 1 A 1 A01 0.659 STD1 STD NA 25 <NA> 0.5545
#> 2 A 2 A02 0.649 STD1 STD NA 25 <NA> 0.5445
#> 3 A 3 A03 0.598 S1 Sample 5 0 Taxol 0.4935
#> 4 A 4 A04 0.601 S1 Sample 5 0 Taxol 0.4965
#> 5 A 5 A05 0.541 S1 Sample 5 0 Taxol 0.4365
#> 6 A 6 A06 0.553 S1 Sample 5 0 Taxol 0.4485
#> estimated
#> 1 26.39687
#> 2 24.01751
#> 3 18.68524
#> 4 18.88785
#> 5 15.70869
#> 6 16.23867
Calculations for standard curve : linear regression curve
fit2<-stats::lm(formula = con ~ OD,data = std)
ggplot2::ggplot(std, ggplot2::aes(x=OD,y=con))+
ggplot2::ggtitle(paste("Std Curve:", nickname))+
ggplot2::geom_point(color="red",size=2)+
ggplot2::geom_line(data = ggplot2::fortify(fit2),ggplot2::aes(x=OD,y=.fitted),
colour="seagreen4",size=1)+
ggplot2::theme_bw()
To evaluate linear fitting
conpred<-estimate(std,colname="OD",fitformula=fit2,method="linear")
compare<-conpred[,c(1,3)]
corAccuracy<-cor(compare)[1,2]
corAccuracy ## Correlation accuracy of regression line
#> [1] 0.9932859
summary(fit2)
#>
#> Call:
#> stats::lm(formula = con ~ OD, data = std)
#>
#> Residuals:
#> 1 2 3 4 5 6 7
#> 0.86249 0.39094 -0.48415 0.06559 -0.16993 -1.92329 1.25834
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -1.1753 0.6078 -1.934 0.111
#> OD 45.3448 2.3618 19.199 7.07e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Residual standard error: 1.135 on 5 degrees of freedom
#> Multiple R-squared: 0.9866, Adjusted R-squared: 0.9839
#> F-statistic: 368.6 on 1 and 5 DF, p-value: 7.069e-06
To estimate the values based on linear curve
estimated_lr<-estimate(data_DF,colname="blankminus",fitformula=fit2,method="linear")
head(estimated_lr)
#> row col position value type id dilution concentration compound blankminus
#> 1 A 1 A01 0.659 STD1 STD NA 25 <NA> 0.5545
#> 2 A 2 A02 0.649 STD1 STD NA 25 <NA> 0.5445
#> 3 A 3 A03 0.598 S1 Sample 5 0 Taxol 0.4935
#> 4 A 4 A04 0.601 S1 Sample 5 0 Taxol 0.4965
#> 5 A 5 A05 0.541 S1 Sample 5 0 Taxol 0.4365
#> 6 A 6 A06 0.553 S1 Sample 5 0 Taxol 0.4485
#> estimated
#> 1 23.96838
#> 2 23.51493
#> 3 21.20234
#> 4 21.33838
#> 5 18.61769
#> 6 19.16183
For multiply estimated by dilution
estimated_lr$estimated2 <- estimated_lr$estimated * estimated_lr$dilution
head(estimated_lr)
#> row col position value type id dilution concentration compound blankminus
#> 1 A 1 A01 0.659 STD1 STD NA 25 <NA> 0.5545
#> 2 A 2 A02 0.649 STD1 STD NA 25 <NA> 0.5445
#> 3 A 3 A03 0.598 S1 Sample 5 0 Taxol 0.4935
#> 4 A 4 A04 0.601 S1 Sample 5 0 Taxol 0.4965
#> 5 A 5 A05 0.541 S1 Sample 5 0 Taxol 0.4365
#> 6 A 6 A06 0.553 S1 Sample 5 0 Taxol 0.4485
#> estimated estimated2
#> 1 23.96838 NA
#> 2 23.51493 NA
#> 3 21.20234 106.01172
#> 4 21.33838 106.69189
#> 5 18.61769 93.08844
#> 6 19.16183 95.80913
For summarising the “estimated_lr” based on “id” and “type”.
result<-dfsummary(estimated_lr,"estimated2",c("id","type"),
c("STD","Blank"),"plate1", rm="FALSE",
param=c(strict="FALSE",cutoff=40,n=12))
#> F1
#> F2
result
#> id type label N Mean SD CV
#> 1 Sample S1 plate1 10 97.804 7.325 7.49
#> 2 Sample S2 plate1 10 92.680 5.703 6.15
#> 3 Sample S3 plate1 10 78.351 10.970 14.00
#> 4 Sample S4 plate1 10 66.811 5.128 7.68
#> 5 Sample S5 plate1 10 60.213 6.797 11.29
#> 6 Sample S6 plate1 10 34.843 5.329 15.30
#> 7 Sample S7 plate1 10 31.850 4.094 12.85
#> 8 Sample S8 plate1 10 38.062 5.869 15.42
For the t test (S1 as control)
pval<-pvalue(result, control="S1", sigval=0.05)
head(pval)
#> id type label N Mean SD CV pvalue significance
#> 1 Sample S1 plate1 10 97.804 7.325 7.49 control
#> 2 Sample S2 plate1 10 92.680 5.703 6.15 0.049 Yes
#> 3 Sample S3 plate1 10 78.351 10.970 14.00 < 0.001 Yes
#> 4 Sample S4 plate1 10 66.811 5.128 7.68 < 0.001 Yes
#> 5 Sample S5 plate1 10 60.213 6.797 11.29 < 0.001 Yes
#> 6 Sample S6 plate1 10 34.843 5.329 15.30 < 0.001 Yes
This example data contain result of dose response of few drugs (drug1, drug2, drug3, drug4) at three concentrations (C1,C2,C3) from two different cell lines (hepg2 and huh7)
Reading the spectrophotometer readings in csv file
rawdata2<-read.csv("384.csv",stringsAsFactors = FALSE,strip.white = TRUE,
na.strings = c("NA",""),header = TRUE,skip=1)
dim(rawdata2)
head(rawdata2)
#> [1] 16 25
#> X X1 X2 X3 X4 X5 X6 X7 X8 X9 X10 X11 X12 X13 X14 X15 X16 X17 X18 X19
#> 1 A 0.1 0.1 0.2 0.3 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.1 0.2 0.3 0.4 0.5 0.6 0.5 0.6
#> 2 B 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.3
#> 3 C 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.4
#> 4 D 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.5
#> 5 E 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.6
#> 6 F 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.7
#> X20 X21 X22 X23 X24
#> 1 0.7 0.8 0.9 1.0 0.1
#> 2 0.3 0.3 0.3 0.3 0.3
#> 3 0.4 0.4 0.4 0.4 0.4
#> 4 0.5 0.5 0.5 0.5 0.5
#> 5 0.2 0.2 0.3 0.3 0.2
#> 6 0.3 0.3 0.4 0.4 0.3
Reading metadata
metadata2<-read.csv("metafile_384_plate.csv",stringsAsFactors = FALSE,strip.white = TRUE,
na.strings = c("NA",""),header = TRUE)
head(metadata2)
#> row col position cell compound concentration type dilution
#> 1 A 1 A01 hepg2 drug1 B blank1 NA
#> 2 A 2 A02 hepg2 drug1 C1 treated 10
#> 3 A 3 A03 hepg2 drug1 C1 treated 10
#> 4 A 4 A04 hepg2 drug1 C1 treated 10
#> 5 A 5 A05 hepg2 drug1 C1 treated 10
#> 6 A 6 A06 hepg2 drug1 C1 treated 10
For appropriate naming of rows and column of rawdata2
rawdata2<-data2plateformat(rawdata2,platetype = 384)
head(rawdata2)
#> 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
#> A 0.1 0.1 0.2 0.3 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.1 0.2 0.3 0.4 0.5 0.6 0.5 0.6
#> B 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.3
#> C 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.4
#> D 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.5
#> E 0.5 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.5 0.5 0.5 0.6 0.6 0.6 0.6 0.6
#> F 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.7
#> 20 21 22 23 24
#> A 0.7 0.8 0.9 1.0 0.1
#> B 0.3 0.3 0.3 0.3 0.3
#> C 0.4 0.4 0.4 0.4 0.4
#> D 0.5 0.5 0.5 0.5 0.5
#> E 0.2 0.2 0.3 0.3 0.2
#> F 0.3 0.3 0.4 0.4 0.3
For heat map
data_DF2<- dplyr::inner_join(OD_df2,metadata2,by=c("row","col","position"))
head(data_DF2)
#> row col position value cell compound concentration type dilution
#> 1 A 1 A01 0.1 hepg2 drug1 B blank1 NA
#> 2 A 2 A02 0.1 hepg2 drug1 C1 treated 10
#> 3 A 3 A03 0.2 hepg2 drug1 C1 treated 10
#> 4 A 4 A04 0.3 hepg2 drug1 C1 treated 10
#> 5 A 5 A05 0.2 hepg2 drug1 C1 treated 10
#> 6 A 6 A06 0.2 hepg2 drug1 C1 treated 10
For categorical view of cells
For categorical view of compound
The data contain separate blanks for drug1, drug2, drug3 and drug4 (blank1, blank2,blank3 and blank 4 respectively)
data_blk<-reduceblank(data_DF2,
x_vector=c("drug1","drug2","drug3","drug4"),
blank_vector = c("blank1","blank2","blank3","blank4"), "value")
dim(data_blk)
#> [1] 384 10
head(data_blk)
#> row col position value cell compound concentration type dilution
#> 1 A 1 A01 0.1 hepg2 drug1 B blank1 NA
#> 2 A 2 A02 0.1 hepg2 drug1 C1 treated 10
#> 3 A 3 A03 0.2 hepg2 drug1 C1 treated 10
#> 4 A 4 A04 0.3 hepg2 drug1 C1 treated 10
#> 5 A 5 A05 0.2 hepg2 drug1 C1 treated 10
#> 6 A 6 A06 0.2 hepg2 drug1 C1 treated 10
#> blankminus
#> 1 -0.4
#> 2 -0.4
#> 3 -0.3
#> 4 -0.2
#> 5 -0.3
#> 6 -0.3
For summarising the result in the order cell, compound,concentration,type and by omitting blanks.
result2<-dfsummary(data_blk,"blankminus",
c("cell","compound","concentration","type"),
c("blank1","blank2","blank3","blank4"),
nickname="384well",
rm="FALSE",param=c(strict="FALSE",cutoff=40,n=12))
#> F1
#> F2
#> F3
#> F4
head (result2)
#> cell compound concentration type label N Mean SD CV
#> 1 hepg2 drug1 C1 treated 384well 7 -0.300 0.058 19.25
#> 2 hepg2 drug1 C1 untreated 384well 8 -0.250 0.053 21.38
#> 3 hepg2 drug1 C2 treated 384well 8 -0.238 0.130 54.84
#> 4 hepg2 drug1 C2 untreated 384well 8 -0.238 0.052 21.79
#> 5 hepg2 drug1 C3 treated 384well 8 0.150 0.278 185.16
#> 6 hepg2 drug1 C3 untreated 384well 8 -0.200 0.000 0.00
dim (result2)
#> [1] 48 9
For summarising the result in the order cell, compound and concentration and by omitting blanks (all blanks are marked as “B” in concentration), drug 2 and huh7.
result3<-dfsummary(data_blk,"blankminus",
c("cell","compound","concentration"),
c("B","drug2","huh7"),
nickname="",
rm="FALSE",param=c(strict="FALSE",cutoff=40,n=12))
#> F1
#> F2
#> F3
head (result3)
#> cell compound concentration label N Mean SD CV
#> 1 hepg2 drug1 C1 15 -0.273 0.059 21.72
#> 2 hepg2 drug1 C2 16 -0.238 0.096 40.31
#> 3 hepg2 drug1 C3 16 -0.025 0.262 1048.17
#> 4 hepg2 drug3 C1 15 0.207 0.070 34.05
#> 5 hepg2 drug3 C2 16 0.212 0.072 33.83
#> 6 hepg2 drug3 C3 16 0.025 0.191 765.94
dim (result3)
#> [1] 9 8
For t-test
pvalue<-pvalue(result3,"C3",sigval=0.05)
pvalue
#> cell compound concentration label N Mean SD CV pvalue
#> 1 hepg2 drug1 C3 16 -0.025 0.262 1048.17 control
#> 2 hepg2 drug1 C1 15 -0.273 0.059 21.72 < 0.001
#> 3 hepg2 drug1 C2 16 -0.238 0.096 40.31 0.0024
#> 4 hepg2 drug3 C3 16 0.025 0.191 765.94 control
#> 5 hepg2 drug3 C1 15 0.207 0.070 34.05 < 0.001
#> 6 hepg2 drug3 C2 16 0.212 0.072 33.83 < 0.001
#> 7 hepg2 drug4 C3 16 0.063 0.266 424.83 control
#> 8 hepg2 drug4 C1 15 0.207 0.070 34.05 0.0258
#> 9 hepg2 drug4 C2 16 0.213 0.072 33.83 0.0187
#> significance
#> 1
#> 2 Yes
#> 3 Yes
#> 4
#> 5 Yes
#> 6 Yes
#> 7
#> 8 Yes
#> 9 Yes