huang 0 like s 620 view s
Tags: Tools
# ------------- Agilent microarray analysis -------------
# Load the necessary packages
library(Biobase)
library(limma)
library(agilp)
library(ggplot2)
# Read in the microarray data
targets <- readTargets("Metadata_Modified.txt", row.names="sample")
rawData <- read.maimages(targets$FileName,source = "agilent.median",green.only = TRUE)
# Pre-process the data using quantile normalization
normalizedData <- normalizeBetweenArrays(rawData, method="quantile")
## Perform PCA on the data
#pca <- prcomp(t(exprs(normalizedData)))
#pca_df <- as.data.frame(pca$x)
#pca_df$Condition <- sampleInfo$condition
#
## Create a scatter plot of the PCA results
#ggplot(pca_df, aes(x=PC1, y=PC2, color=Condition)) +
# geom_point(size=3) +
# xlab(paste0("PC1 (", round(summary(pca)$importance[2,1]*100), "%)")) +
# ylab(paste0("PC2 (", round(summary(pca)$importance[2,2]*100), "%)")) +
# ggtitle("PCA Plot") +
# theme(plot.title = element_text(hjust = 0.5))
# Extract the sample information and create a design matrix
sampleInfo <- read.table("Metadata_Modified.txt", header=TRUE)
design <- model.matrix(~ 0 + sampleInfo$condition)
# Fit a linear model to the data
fit <- lmFit(normalizedData, design)
# Compute moderated t-statistics and p-values
fit <- eBayes(fit)
# Identify differentially expressed genes using a specified FDR threshold
diffGenes <- topTable(fit, coef=1, adjust.method="fdr", sort.by="B", n=Inf)
# Output the results
write.csv(diffGenes, file="diffGenes.csv")
# Create a heatmap of the differentially expressed genes
heatmapExprs <- exprs(normalizedData)[diffGenes$ID,]
heatmapCol <- colorRampPalette(c("blue","white","red"))(50)
heatmap(heatmapExprs, Colv=NA, col=heatmapCol, scale="row", margins=c(10,10))
# ------------- Affymetrix microarray analysis -------------
# Load the necessary packages
library(limma)
library(affy)
# Read in the microarray data
rawData <- ReadAffy()
# Pre-process the data using RMA normalization
normalizedData <- expresso(rawData, normalize.method="rma")
# Extract the sample information and create a design matrix
sampleInfo <- read.table("sample_info.txt", header=TRUE)
design <- model.matrix(~ 0 + sampleInfo$Group)
# Fit a linear model to the data
fit <- lmFit(normalizedData, design)
# Compute moderated t-statistics and p-values
fit <- eBayes(fit)
# Identify differentially expressed genes using a specified FDR threshold
diffGenes <- topTable(fit, coef=1, adjust.method="fdr", sort.by="B", n=Inf)
# Perform functional enrichment analysis on the differentially expressed genes
enrichment <- enrichKEGG(diffGenes$ID, organism="hsa", pvalueCutoff=0.05)
# Print the results
print(diffGenes)
print(enrichment)
点赞本文的读者
还没有人对此文章表态
没有评论
Run GAMOLA2 on MT880870-MT880872
Why are all adjusted p-values are the same in the DESeq2 result?
Bactopia: a Flexible Pipeline for Complete Analysis of Bacterial Genomes
© 2023 XGenes.com Impressum