Import file

x <- read_excel("C://HB-histology-July-2021/measurements.xlsx", 
                col_types = c("text", "text", "numeric", 
                              "numeric", "numeric", "numeric",
                              "numeric", "numeric", "numeric"))

DATA MANIPULATION

Average measures by sample*ep_cell_depth

total <- x %>%                              # Specify data frame
  group_by(sample, genotype, ep_cell_depth) %>%           # Specify group indicator
  summarise_at(vars(total_depth),                         # Specify column
               list(total_mean = mean)) 

ep <- x %>%                              # Specify data frame
  group_by(sample, genotype, ep_cell_depth) %>%           # Specify group indicator
  summarise_at(vars(ep_depth),                         # Specify column
               list(ep_mean = mean)) 

mesen <- x %>%                              # Specify data frame
  group_by(sample, genotype, ep_cell_depth) %>%           # Specify group indicator
  summarise_at(vars(mesen_depth),                         # Specify column
               list(mesen_mean = mean)) 

dense <- x %>%                              # Specify data frame
  group_by(sample, genotype, ep_cell_depth) %>%           # Specify group indicator
  summarise_at(vars(dense_depth),                         # Specify column
               list(dense_mean = mean)) 

y <- left_join(total, ep, by = c("sample", "genotype", "ep_cell_depth"))
y <- left_join(y, mesen, by = c("sample", "genotype", "ep_cell_depth"))
y <- left_join(y, dense, by = c("sample", "genotype", "ep_cell_depth"))


y <- within(y, ep_ratio <- ep_mean / total_mean)
y <- within(y, mesen_ratio <- mesen_mean / total_mean)
y <- within(y, dens_ratio <- dense_mean / total_mean)

#split y by epithelium cell depth and bind to form innerHB, outerHB and polledHB+FS

split<- split(y, y$ep_cell_depth)

outer <- bind_rows(split$"1",split$"1.5", split$"2")
outer<- split(outer, outer$genotype)
outerpp<-bind_rows(outer$"pp")
PP<-bind_rows(outer$"PP")

inner <- bind_rows(split$"7", split$"7.5", split$"8", split$"8.5", split$"9", split$"9.5", split$"10")


outerpp$i <-"OuterHB"
inner$i <-"InnerHB"
PP$i <- "PolledHB+FS"

four <-bind_rows(inner,outerpp,PP)

TEST FOR NORMALITY (MEANS)

OuterHB

total

ggdensity(outerpp$total_mean)

ggqqplot(outerpp$total_mean)

shapiro.test(outerpp$total_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$total_mean
## W = 0.92967, p-value = 0.3019

epithelium mean

ggdensity(outerpp$ep_mean)

ggqqplot(outerpp$ep_mean)

shapiro.test(outerpp$ep_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$ep_mean
## W = 0.90613, p-value = 0.1385

mesenchyme mean

ggdensity(outerpp$mesen_mean)

ggqqplot(outerpp$mesen_mean)

shapiro.test(outerpp$mesen_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$mesen_mean
## W = 0.93559, p-value = 0.3648

condense cell mean

ggdensity(outerpp$dense_mean)

ggqqplot(outerpp$dense_mean)

shapiro.test(outerpp$dense_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$dense_mean
## W = 0.95016, p-value = 0.5632

inner

total mean

ggdensity(inner$total_mean)

ggqqplot(inner$total_mean)

shapiro.test(inner$total_mean)#normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$total_mean
## W = 0.9387, p-value = 0.4403

epithelium mean

ggdensity(inner$ep_mean)

ggqqplot(inner$ep_mean)

shapiro.test(inner$ep_mean) #not normal, outlier

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$ep_mean
## W = 0.83541, p-value = 0.01853

mesenchyme mean

ggdensity(inner$mesen_mean)

ggqqplot(inner$mesen_mean)

shapiro.test(inner$mesen_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$mesen_mean
## W = 0.97789, p-value = 0.9677

condense cell mean

ggdensity(inner$dense_mean)

ggqqplot(inner$dense_mean)

shapiro.test(inner$dense_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$dense_mean
## W = 0.95529, p-value = 0.6801

PP

total mean

ggdensity(PP$total_mean)

ggqqplot(PP$total_mean)

shapiro.test(PP$total_mean)#normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$total_mean
## W = 0.90174, p-value = 0.3843

epithelium mean

ggdensity(PP$ep_mean)

ggqqplot(PP$ep_mean)

shapiro.test(PP$ep_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$ep_mean
## W = 0.92839, p-value = 0.5677

mesenchyme mean

ggdensity(PP$mesen_mean)

ggqqplot(PP$mesen_mean)

shapiro.test(PP$mesen_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$mesen_mean
## W = 0.88679, p-value = 0.3017

condense cell mean

ggdensity(PP$dense_mean)

ggqqplot(PP$dense_mean)

shapiro.test(PP$dense_mean) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$dense_mean
## W = 0.97501, p-value = 0.9242

Epithelium mean is not normally distributed, will test normality when converted to a proportion of the total depth.

TEST FOR NORMALITY (RATIOS)

OuterHB

epithelium

ggdensity(outerpp$ep_ratio)

ggqqplot(outerpp$ep_ratio)

shapiro.test(outerpp$ep_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$ep_ratio
## W = 0.90847, p-value = 0.1497

mesenchyme

ggdensity(outerpp$mesen_ratio)

ggqqplot(outerpp$mesen_ratio)

shapiro.test(outerpp$mesen_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$mesen_ratio
## W = 0.9607, p-value = 0.7346

###condesnsed cells

ggdensity(outerpp$dens_ratio)

ggqqplot(outerpp$dens_ratio)

shapiro.test(outerpp$dens_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  outerpp$dens_ratio
## W = 0.98922, p-value = 0.9992

inner

epithelium

ggdensity(inner$ep_ratio)

ggqqplot(inner$ep_ratio)

shapiro.test(inner$ep_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$ep_ratio
## W = 0.97074, p-value = 0.903

mesenchyme

ggdensity(inner$mesen_ratio)

ggqqplot(inner$mesen_ratio)

shapiro.test(inner$mesen_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$mesen_ratio
## W = 0.94951, p-value = 0.5908

condesnsed cells

ggdensity(inner$dens_ratio)

ggqqplot(inner$dens_ratio)

shapiro.test(inner$dens_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  inner$dens_ratio
## W = 0.92757, p-value = 0.3167

PP

epithelium

ggdensity(PP$ep_ratio)

ggqqplot(PP$ep_ratio)

shapiro.test(PP$ep_ratio) #not normal, outlier

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$ep_ratio
## W = 0.63467, p-value = 0.001176

mesenchyme

ggdensity(PP$mesen_ratio)

ggqqplot(PP$mesen_ratio)

shapiro.test(PP$mesen_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$mesen_ratio
## W = 0.83312, p-value = 0.1142

condesnsed cells

ggdensity(PP$dens_ratio)

ggqqplot(PP$dens_ratio)

shapiro.test(PP$dens_ratio) #normal

## 
##  Shapiro-Wilk normality test
## 
## data:  PP$dens_ratio
## W = 0.89082, p-value = 0.3225

Therefore the data is not normal and a non-parametric test must be used.

VISUALISATION

boxplots

#SET STATS COMPARISIONS

#t <- list(c('InnerHB', 'OuterHB'), c('InnerHB', 'PolledHB+FS'), c('OuterHB', 'PolledHB+FS'))



g1 <- four%>%
  ggplot( aes(x=i, y=total_mean, fill = genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
                    # label = "p.signif", 
                    # symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
                                      #  symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Total depth") + 
  xlab("Position") +
  ylab("Depth (um)") +
  theme_light() +
  theme(plot.title = element_text(hjust = 0.5))

g2 <- four%>%
  ggplot( aes(x=i, y=ep_ratio, fill = genotype))+
  geom_boxplot(width = 0.5)+
 # stat_compare_means(comparison = t,
                     #label = "p.signif", 
                     #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
                                        #symbols = c("****", "***", "**", "*", "ns")))+
  #ggtitle("Epithelium depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_light()  +
  theme(plot.title = element_text(hjust = 0.5))

g3 <- four%>%
  ggplot( aes(x=i, y=mesen_ratio, fill = genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
                     #label = "p.signif", 
                     #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
                                        #symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Mesenchyme depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_light()  +
  theme(plot.title = element_text(hjust = 0.5))

g4 <- four%>%
  ggplot( aes(x=i, y=dens_ratio, fill = genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
                     #label = "p.signif", 
                     #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
                                          #symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Condense cell depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_light()  +
  theme(plot.title = element_text(hjust = 0.5))


g<-grid.arrange(g1, g2, g3, g4, nrow=2)

DESCRIPTIVE STATS

total

ds_total <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(total_mean, na.rm = TRUE),
    median = median(total_mean, na.rm = TRUE),
    sd = sd(total_mean, na.rm = TRUE)
  )
ds_total

## # A tibble: 3 x 5
##   i           count  mean median    sd
##   <chr>       <int> <dbl>  <dbl> <dbl>
## 1 InnerHB        13  225.   228.  42.5
## 2 OuterHB        14  146.   151.  22.2
## 3 PolledHB+FS     6  155.   163.  41.5

write.csv(ds_total,"C:/HB-histology-July-2021/ds_total.csv", row.names = FALSE)

epithlium

ds_ep <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(ep_mean, na.rm = TRUE),
    median = median(ep_mean, na.rm = TRUE),
    sd = sd(ep_mean, na.rm = TRUE)
  )
ds_ep

## # A tibble: 3 x 5
##   i           count  mean median    sd
##   <chr>       <int> <dbl>  <dbl> <dbl>
## 1 InnerHB        13 89.0   83.3  19.5 
## 2 OuterHB        14 13.0   11.9   4.45
## 3 PolledHB+FS     6  9.78   9.03  3.08

write.csv(ds_ep,"C:/HB-histology-July-2021/ds_ep.csv", row.names = FALSE)

mesenchyme

ds_mesen <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(mesen_mean, na.rm = TRUE),
    median = median(mesen_mean, na.rm = TRUE),
    sd = sd(mesen_mean, na.rm = TRUE)
  )
ds_mesen

## # A tibble: 3 x 5
##   i           count  mean median    sd
##   <chr>       <int> <dbl>  <dbl> <dbl>
## 1 InnerHB        13  70.0   71.2  20.1
## 2 OuterHB        14 117.   114.   20.6
## 3 PolledHB+FS     6 132.   142.   47.4

write.csv(ds_mesen,"C:/HB-histology-July-2021/ds_mesen.csv", row.names = FALSE)

condense

ds_dense <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(dense_mean, na.rm = TRUE),
    median = median(dense_mean, na.rm = TRUE),
    sd = sd(dense_mean, na.rm = TRUE)
  )
ds_dense

## # A tibble: 3 x 5
##   i           count  mean median    sd
##   <chr>       <int> <dbl>  <dbl> <dbl>
## 1 InnerHB        13  65.6   65.3 16.9 
## 2 OuterHB        14  15.5   14.5  9.77
## 3 PolledHB+FS     6  13.0   12.7 10.1

write.csv(ds_dense,"C:/HB-histology-July-2021/ds_dense.csv", row.names = FALSE)

Mean ratio (epithelium, mesenchyme and condensed cell layer) calculated for Inner, Outer, PolledHB+FS

inner

i<- split(four, four$i)

inner <- i$"InnerHB"
head(inner)

## # A tibble: 6 x 11
## # Groups:   sample, genotype [6]
##   sample genotype ep_cell_depth total_mean ep_mean mesen_mean dense_mean
##   <chr>  <chr>            <dbl>      <dbl>   <dbl>      <dbl>      <dbl>
## 1 532HB  pp                   7       165.    76.7       50.1       37.9
## 2 546HB  pp                   7       228.    83.3       76.1       65.3
## 3 581HB  pp                   7       280.    80.0      112.        83.8
## 4 618HB  pp                   7       175.    71.1       60.3       43.7
## 5 668HB  pp                   7       211.    79.7       53.2       76.9
## 6 736HB  pp                   7       246.    85.0       83.8       77.3
## # ... with 4 more variables: ep_ratio <dbl>, mesen_ratio <dbl>,
## #   dens_ratio <dbl>, i <chr>

inner_stats <- group_by(inner, sample) %>%
  summarise(count = n(),
            mean_total_depth = mean(total_mean, na.rm = TRUE),
    mean_ep_ratio = mean(ep_ratio, na.rm = TRUE),
    mean_mesen_ratio = mean(mesen_ratio, na.rm = TRUE),
    mean_dens_ratio = mean(dens_ratio, na.rm = TRUE)
  )
head(inner_stats)

## # A tibble: 6 x 6
##   sample count mean_total_depth mean_ep_ratio mean_mesen_ratio mean_dens_ratio
##   <chr>  <int>            <dbl>         <dbl>            <dbl>           <dbl>
## 1 532HB      2             156.         0.460            0.277           0.268
## 2 546HB      2             230.         0.367            0.342           0.280
## 3 581HB      1             280.         0.286            0.399           0.299
## 4 618HB      4             230.         0.454            0.301           0.246
## 5 668HB      2             213.         0.380            0.244           0.382
## 6 736HB      2             262.         0.343            0.338           0.321

outer

outer<- i$"OuterHB"
outer_stats <- group_by(outer, sample) %>%
  summarise(count = n(),
            mean_total_depth = mean(total_mean, na.rm = TRUE),
            mean_ep_ratio = mean(ep_ratio, na.rm = TRUE),
            mean_mesen_ratio = mean(mesen_ratio, na.rm = TRUE),
            mean_dens_ratio = mean(dens_ratio, na.rm = TRUE)
  )
head(outer_stats)

## # A tibble: 6 x 6
##   sample count mean_total_depth mean_ep_ratio mean_mesen_ratio mean_dens_ratio
##   <chr>  <int>            <dbl>         <dbl>            <dbl>           <dbl>
## 1 532HB      3             137.        0.0880            0.754          0.141 
## 2 546HB      2             169.        0.0970            0.758          0.145 
## 3 581HB      3             159.        0.0639            0.891          0.0421
## 4 618FS      2             105.        0.0761            0.828          0.0906
## 5 618HB      2             140.        0.0957            0.740          0.151 
## 6 668HB      2             163.        0.122             0.798          0.0799

polledHB+FS

polled<-i$"PolledHB+FS"
polled_stats <- group_by(polled, sample) %>%
  summarise(count = n(),
            mean_total_depth = mean(total_mean, na.rm = TRUE),
            mean_ep_ratio = mean(ep_ratio, na.rm = TRUE),
            mean_mesen_ratio = mean(mesen_ratio, na.rm = TRUE),
            mean_dens_ratio = mean(dens_ratio, na.rm = TRUE)
  )

head(polled_stats)

## # A tibble: 3 x 6
##   sample count mean_total_depth mean_ep_ratio mean_mesen_ratio mean_dens_ratio
##   <chr>  <int>            <dbl>         <dbl>            <dbl>           <dbl>
## 1 667HB      2             108.        0.0726            0.710          0.218 
## 2 694FS      2             193.        0.0696            0.917          0.0125
## 3 709HB      2             163.        0.0515            0.869          0.0761

Combine data

all_mean_data <- rbind(inner_stats, outer_stats, polled_stats)
all_mean_data$i <- c("innerHB","innerHB","innerHB","innerHB","innerHB","innerHB",
                     "outerHB","outerHB", "outerHB", "outerHB","outerHB","outerHB",
                     "polledHB+FS","polledHB+FS","polledHB+FS")

all_mean_data$genotype <- c("horned","horned","horned","horned","horned","horned",
                     "horned","horned", "horned", "horned","horned","horned",
                     "polled","polled","polled")

NORMALITY TEST OF MEAN DATA

all data

total depth

ggdensity(all_mean_data$mean_total_depth) #normal

ggqqplot(all_mean_data$mean_total_depth)

shapiro.test(all_mean_data$mean_total_depth)

## 
##  Shapiro-Wilk normality test
## 
## data:  all_mean_data$mean_total_depth
## W = 0.94759, p-value = 0.4874

epithelium ratio

ggdensity(all_mean_data$mean_ep_ratio) #not normal

ggqqplot(all_mean_data$mean_ep_ratio)

shapiro.test(all_mean_data$mean_ep_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  all_mean_data$mean_ep_ratio
## W = 0.79413, p-value = 0.003094

mesechyme ratio

ggdensity(all_mean_data$mean_mesen_ratio) #not normal

ggqqplot(all_mean_data$mean_mesen_ratio)

shapiro.test(all_mean_data$mean_mesen_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  all_mean_data$mean_mesen_ratio
## W = 0.83628, p-value = 0.01117

condensed cell ratio

ggdensity(all_mean_data$mean_dens_ratio) #normal

ggqqplot(all_mean_data$mean_dens_ratio)

shapiro.test(all_mean_data$mean_dens_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  all_mean_data$mean_dens_ratio
## W = 0.95572, p-value = 0.6186

inner

total depth

ggdensity(inner_stats$mean_total_depth) #normal

ggqqplot(inner_stats$mean_total_depth)

shapiro.test(inner_stats$mean_total_depth)

## 
##  Shapiro-Wilk normality test
## 
## data:  inner_stats$mean_total_depth
## W = 0.94263, p-value = 0.6804

epithelium ratio

ggdensity(inner_stats$mean_ep_ratio) #normal

ggqqplot(inner_stats$mean_ep_ratio)

shapiro.test(inner_stats$mean_ep_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  inner_stats$mean_ep_ratio
## W = 0.93439, p-value = 0.6144

mesechyme ratio

ggdensity(inner_stats$mean_mesen_ratio) #normal

ggqqplot(inner_stats$mean_mesen_ratio)

shapiro.test(inner_stats$mean_mesen_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  inner_stats$mean_mesen_ratio
## W = 0.98228, p-value = 0.9623

condensed cell ratio

ggdensity(inner_stats$mean_dens_ratio) #normal

ggqqplot(inner_stats$mean_dens_ratio)

shapiro.test(inner_stats$mean_dens_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  inner_stats$mean_dens_ratio
## W = 0.94005, p-value = 0.6596

outer

total depth

ggdensity(outer_stats$mean_total_depth) #normal

ggqqplot(outer_stats$mean_total_depth)

shapiro.test(outer_stats$mean_total_depth)

## 
##  Shapiro-Wilk normality test
## 
## data:  outer_stats$mean_total_depth
## W = 0.90172, p-value = 0.3842

epithelium ratio

ggdensity(outer_stats$mean_ep_ratio) #normal

ggqqplot(outer_stats$mean_ep_ratio)

shapiro.test(outer_stats$mean_ep_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  outer_stats$mean_ep_ratio
## W = 0.97226, p-value = 0.9072

mesechyme ratio

ggdensity(outer_stats$mean_mesen_ratio) #normal

ggqqplot(outer_stats$mean_mesen_ratio)

shapiro.test(outer_stats$mean_mesen_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  outer_stats$mean_mesen_ratio
## W = 0.89708, p-value = 0.3569

condensed cell ratio

ggdensity(outer_stats$mean_dens_ratio) #normal

ggqqplot(outer_stats$mean_dens_ratio)

shapiro.test(outer_stats$mean_dens_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  outer_stats$mean_dens_ratio
## W = 0.87897, p-value = 0.2644

polled data

total depth

ggdensity(polled_stats$mean_total_depth) #normal

ggqqplot(polled_stats$mean_total_depth)

shapiro.test(polled_stats$mean_total_depth)

## 
##  Shapiro-Wilk normality test
## 
## data:  polled_stats$mean_total_depth
## W = 0.97265, p-value = 0.6827

epithelium ratio

ggdensity(polled_stats$mean_ep_ratio) #normal

ggqqplot(polled_stats$mean_ep_ratio)

shapiro.test(polled_stats$mean_ep_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  polled_stats$mean_ep_ratio
## W = 0.85398, p-value = 0.2512

mesechyme ratio

ggdensity(polled_stats$mean_mesen_ratio) # normal

ggqqplot(polled_stats$mean_mesen_ratio)

shapiro.test(polled_stats$mean_mesen_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  polled_stats$mean_mesen_ratio
## W = 0.91274, p-value = 0.4273

condensed cell ratio

ggdensity(polled_stats$mean_dens_ratio) #normal

ggqqplot(polled_stats$mean_dens_ratio)

shapiro.test(polled_stats$mean_dens_ratio)

## 
##  Shapiro-Wilk normality test
## 
## data:  polled_stats$mean_dens_ratio
## W = 0.95377, p-value = 0.5861

BOXPLOTS OF AVERAGED RATIOS

all_mean_data %>% head

## # A tibble: 6 x 8
##   sample count mean_total_depth mean_ep_ratio mean_mesen_ratio mean_dens_ratio
##   <chr>  <int>            <dbl>         <dbl>            <dbl>           <dbl>
## 1 532HB      2             156.         0.460            0.277           0.268
## 2 546HB      2             230.         0.367            0.342           0.280
## 3 581HB      1             280.         0.286            0.399           0.299
## 4 618HB      4             230.         0.454            0.301           0.246
## 5 668HB      2             213.         0.380            0.244           0.382
## 6 736HB      2             262.         0.343            0.338           0.321
## # ... with 2 more variables: i <chr>, genotype <chr>

g1 <- all_mean_data%>%
  ggplot( aes(x=i, y=mean_total_depth, fill=genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
  # label = "p.signif", 
  # symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
  #  symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Total depth") + 
  xlab("Position") +
  ylab("Depth (um)") +
  theme_classic() +
  theme(plot.title = element_text(hjust = 0.5))

g2 <- all_mean_data%>%
  ggplot( aes(x=i, y=mean_ep_ratio, fill=genotype))+
  geom_boxplot(width = 0.5)+
  # stat_compare_means(comparison = t,
  #label = "p.signif", 
  #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
  #symbols = c("****", "***", "**", "*", "ns")))+
  #ggtitle("Epithelium depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_classic()  +
  theme(plot.title = element_text(hjust = 0.5))

g3 <- all_mean_data%>%
  ggplot( aes(x=i, y=mean_mesen_ratio, fill=genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
  #label = "p.signif", 
  #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
  #symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Mesenchyme depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_classic()  +
  theme(plot.title = element_text(hjust = 0.5))

g4 <- all_mean_data%>%
  ggplot( aes(x=i, y=mean_dens_ratio, fill=genotype))+
  geom_boxplot(width = 0.5)+
  #stat_compare_means(comparison = t, 
  #label = "p.signif", 
  #symnum.args = list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), 
  #symbols = c("****", "***", "**", "*", "ns"))) +
  #ggtitle("Condense cell depth") +
  xlab("Position") +
  ylab("Proportion of Total Depth") +
  theme_classic()  +
  theme(plot.title = element_text(hjust = 0.5))


g<-grid.arrange(g1, g2, g3, g4, nrow=2)

STATISTICAL TESTS: UNPAIRED & PAIRED WILCOXON RANK SUM TEST (MANN-WHITNEY U TEST)

comparisons: (‘InnerHB’, ‘OuterHB’) (‘InnerHB’, ‘PolledHB+FS’) (‘OuterHB’, ‘PolledHB+FS’)

total depth

total_mean_innervsouter <- wilcox.test(inner_stats$mean_total_depth, outer_stats$mean_total_depth,
            paired = TRUE, alternative = "two.sided")
total_mean_innervsouter$p.value

## [1] 0.03125

total_mean_innervspolled <- wilcox.test(inner_stats$mean_total_depth, polled_stats$mean_total_depth,
                                        paired = FALSE, alternative = "two.sided")
total_mean_innervspolled$p.value

## [1] 0.0952381

total_mean_outervspolled <- wilcox.test(outer_stats$mean_total_depth, polled_stats$mean_total_depth,
                                         paired = FALSE, alternative = "two.sided")
total_mean_outervspolled$p.value

## [1] 0.7142857

epithelium ratio

ep_ratio_innervsouter <- wilcox.test(inner_stats$mean_ep_ratio, outer_stats$mean_ep_ratio,
                                        paired = TRUE, alternative = "two.sided")
ep_ratio_innervsouter$p.value

## [1] 0.03125

ep_ratio_innervspolled <- wilcox.test(inner_stats$mean_ep_ratio, polled_stats$mean_ep_ratio,
                                         paired = FALSE, alternative = "two.sided")
ep_ratio_innervspolled$p.value

## [1] 0.02380952

ep_ratio_outervspolled <- wilcox.test(outer_stats$mean_ep_ratio, polled_stats$mean_ep_ratio,
                                         paired = FALSE, alternative = "two.sided")
ep_ratio_outervspolled$p.value

## [1] 0.0952381

mesenchyme ratio

mesen_ratio_innervsouter <- wilcox.test(inner_stats$mean_mesen_ratio, outer_stats$mean_mesen_ratio,
                                     paired = TRUE, alternative = "two.sided")
mesen_ratio_innervsouter$p.value

## [1] 0.03125

mesen_ratio_innervspolled <- wilcox.test(inner_stats$mean_mesen_ratio, polled_stats$mean_mesen_ratio,
                                      paired = FALSE, alternative = "two.sided")
mesen_ratio_innervspolled$p.value

## [1] 0.02380952

mesen_ratio_outervspolled <- wilcox.test(outer_stats$mean_mesen_ratio, polled_stats$mean_mesen_ratio,
                                      paired = FALSE, alternative = "two.sided")
mesen_ratio_outervspolled$p.value

## [1] 0.7142857

condensed cells ratio

dens_ratio_innervsouter <- wilcox.test(inner_stats$mean_dens_ratio, outer_stats$mean_dens_ratio,
                                        paired = TRUE, alternative = "two.sided")
dens_ratio_innervsouter$p.value

## [1] 0.03125

dens_ratio_innervspolled <- wilcox.test(inner_stats$mean_dens_ratio, polled_stats$mean_dens_ratio,
                                         paired = FALSE, alternative = "two.sided")
dens_ratio_innervspolled$p.value

## [1] 0.02380952

dens_ratio_outervspolled <- wilcox.test(outer_stats$mean_dens_ratio, polled_stats$mean_dens_ratio,
                                         paired = FALSE, alternative = "two.sided")
dens_ratio_outervspolled$p.value

## [1] 0.7142857

CORRELATIONS

#changing genotype 'value' to have clearer label
z <- within(as.data.frame(x), 
            { genotype[genotype == "PP"] <- "Polled" })
z <- within(as.data.frame(z), 
            { genotype[genotype == "pp"] <- "Horned" })

ggplot(z, aes(x = `total_depth`, y = `ep_depth` )) +  #plots graph with x-axis as weight and y-axis as height
  geom_point(aes(colour = `ep_cell_depth`), size = 4) + # 
  geom_smooth(method = "lm", formula = (y ~ x), se = FALSE) +# add line of best fit
  stat_cor(method = "pearson") +
  xlab("Total Depth") +
  ylab("Epithelium Depth") +
  facet_wrap(~genotype) + #split plot by genotype
  theme_grey()# r value

## Warning: Removed 2 rows containing non-finite values (stat_smooth).

## Warning: Removed 2 rows containing non-finite values (stat_cor).

## Warning: Removed 2 rows containing missing values (geom_point).

COLUMN GRAPHS

calculate mean ratios, sd and se for inner, outer and polledhb+fs

ds_total = column graph for total depth
col_two = column graph for epithelium ratio
col_three = column graph for mesenchyme ratio
col_four = column graph for condensed cell ratio

col_two <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(ep_ratio, na.rm = TRUE),
    sd = sd(ep_ratio, na.rm = TRUE)
  )

col_three <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(mesen_ratio, na.rm = TRUE),
    sd = sd(mesen_ratio, na.rm = TRUE)
  )

col_four <- group_by(four, i) %>%
  summarise(
    count = n(),
    mean = mean(dens_ratio, na.rm = TRUE),
    sd = sd(dens_ratio, na.rm = TRUE)
  )


ds_total$n <- c(6, 6, 3) #add column for sample numbers
ds_total$se <- ds_total$sd/sqrt(ds_total$n)
col_two$n <- c(6, 6, 3)
col_two$se <- col_two$sd/sqrt(col_two$n)
col_three$n <- c(6, 6, 3)
col_three$se <- col_three$sd/sqrt(col_three$n)
col_four$n <- c(6, 6, 3)
col_four$se <- col_four$sd/sqrt(col_four$n)

plots

total depth

ggplot(ds_total, aes(x=i, y=mean))+ 
  geom_col(width = 0.7, fill = "steelblue")+
  geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width=.2,
                                                    position=position_dodge(.9)) +
  labs(x = "Sample", y = "Depth (µm)")+
  ylim(0,300)+
  theme_light(base_size = 18)

epithelium proportion of depth

ggplot(col_two, aes(x=i, y=mean))+ 
  geom_col(width = 0.7, fill = "steelblue")+
  geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width=.2,
                position=position_dodge(.9)) +
  labs(x = "Sample", y = "Proportion of Total Depth") +
  ylim(0,0.5)+
  theme_light(base_size = 18)

mesenchyme proportion of depth

ggplot(col_three, aes(x=i, y=mean))+ 
  geom_col(width = 0.7, fill = "steelblue")+
  geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width=.2,
                position=position_dodge(.9)) +
  labs(x = "Sample", y = "Proportion of Total Depth") +
  ylim(0,1)+
  theme_light(base_size = 18)

condensed cells proportion of depth

ggplot(col_four, aes(x=i, y=mean))+ 
  geom_col(width = 0.7, fill = "steelblue")+
  geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width=.2,
                position=position_dodge(.9)) +
  labs(x = "Sample", y = "Proportion of Total Depth") +
  ylim(0,0.5)+
  theme_light(base_size = 18)

Histology_stats

J Aldersey

03/06/2022

Import file

DATA MANIPULATION

Average measures by sample*ep_cell_depth

TEST FOR NORMALITY (MEANS)

OuterHB

total

epithelium mean

mesenchyme mean

condense cell mean

inner

total mean

epithelium mean

mesenchyme mean

condense cell mean

PP

total mean

epithelium mean

mesenchyme mean

condense cell mean

TEST FOR NORMALITY (RATIOS)

OuterHB

epithelium

mesenchyme

inner

epithelium

mesenchyme

condesnsed cells

PP

epithelium

mesenchyme

condesnsed cells

VISUALISATION

boxplots

DESCRIPTIVE STATS

total

epithlium

mesenchyme

condense

Mean ratio (epithelium, mesenchyme and condensed cell layer) calculated for Inner, Outer, PolledHB+FS

inner

outer

polledHB+FS

Combine data

NORMALITY TEST OF MEAN DATA

all data

total depth

epithelium ratio

mesechyme ratio

condensed cell ratio

inner

total depth

epithelium ratio

mesechyme ratio

condensed cell ratio

outer

total depth

epithelium ratio

mesechyme ratio

condensed cell ratio

polled data

total depth

epithelium ratio

mesechyme ratio

condensed cell ratio

BOXPLOTS OF AVERAGED RATIOS

STATISTICAL TESTS: UNPAIRED & PAIRED WILCOXON RANK SUM TEST (MANN-WHITNEY U TEST)

total depth

epithelium ratio

mesenchyme ratio

condensed cells ratio

CORRELATIONS

COLUMN GRAPHS

calculate mean ratios, sd and se for inner, outer and polledhb+fs

plots

total depth

epithelium proportion of depth

mesenchyme proportion of depth