
#if (!requireNamespace("BiocManager", quietly = TRUE))
#    install.packages("BiocManager")
#BiocManager::install("limma", version = "3.12")

 
  #¡¾A¡¿#############################
  ###  ---different expression  ---### 
  ####################################---#-----------------------------------------------------------------
  
  
  library(limma)
  
  filename.IHC.txt <- paste("FTH1.IHC score.spss.txt",sep = "")
  rt=read.table(filename.IHC.txt,sep="\t",header=T,check.names=F)#row.names=F,
 
  library(tidyverse)
  library(ggplot2)
  library(reshape2)
  library(car)
  library(rstatix)
  
  x <- rt[,1]#[grep(">",rtmiRName$Symble),]
  y<- rt[,2]
  
  data <- data.frame(x = x, y = y)
  data2 <- melt(data)
  
  data3 <- lapply(data, function(x) get_summary_stats(data.frame(x)))
  data3
  # $x
  # # A tibble: 1 x 13
  # variable     n    min   max median    q1    q3   iqr   mad  mean    sd    se    ci
  # <chr>    <dbl>  <dbl> <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
  #   1 x          100 -0.272  4.58   1.94  1.39  2.66  1.26 0.974  2.00  1.02 0.102 0.203
  # 
  # $y
  # # A tibble: 1 x 13
  # variable     n   min   max median    q1    q3   iqr   mad  mean    sd    se    ci
  # <chr>    <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
  #   1 x          100 -1.14  3.17  0.927 0.568  1.45 0.878 0.648  1.01 0.796  0.08 0.158
  data3 <- rbind(data3[[1]], data3[[2]])
  data3[1] <- c("x", "y")
  
  ## Shapiro-Wilk normality test
  lapply(data, function(x) shapiro.test(x))
  # $x
  # 
  # Shapiro-Wilk normality test
  # 
  # data:  x
  # W = 0.98836, p-value = 0.535
  # 
  # 
  # $y
  # 
  # Shapiro-Wilk normality test
  # 
  # data:  x
  # W = 0.98532, p-value = 0.3348
  
  
  #p-Value
  
  leveneTest(value~variable, data = data2)
  
  t.test(value~variable, data = data2, var.equal = T)
  
  
  
  summary(aov(value~variable, data = data2))
  
  
  outpdf=paste("Figure 2F.pdf",sep="")
  pdf(outpdf,width=5,height=5)
  
  ggplot(data2, aes(x = variable, y = value, color = variable, fill = variable)) +
    geom_violin(alpha = 0.2) +
    geom_point(position = position_jitter(0.3)) +
    theme_bw()
  dev.off()
  