Entering edit mode
6 months ago
Natali
•
0
Hi! I can't solve this problem. Here there is my code, can someone help me to find and solve the code? I am stuck because I can't normalise my data, as it says "DGEList function not found", but maybe there are other errors too.
Thank you
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("genefu")
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("org.Hs.eg.db")
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("edgeR")
Annotation hub if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("AnnotationHub")
library(genefu)
library(AnnotationDbi)
library(org.Hs.eg.db)
library(dplyr)
library(edgeR)
#Cargar datos clínicos y leerlos
library(genefu)
library(AnnotationDbi)
library(org.Hs.eg.db)
library(dplyr)
library(edgeR)
mat<-read.table("~/Downloads/BRCA_exp_matrix.tsv",header=TRUE,sep="\t",fill=TRUE)
library(readr)
clinical <- read.table("~/Downloads/clinical_info_TCGA-BRCA.tsv", sep = "\t", na.strings = "NAN", stringsAsFactors = FALSE)
#filtrado de muestras
#Seleccionamos las muestras que sean Luminales, que expresan el receptor de estrogenos o sean ER positivas por IHC
er_query<-!is.na(clinical$er_status_by_ihc)
clinical<-clinical[er_query,]
sum(clinical$er_status_by_ihc=="Positive")
clinical<-clinical[clinical$er_status_by_ihc=="Positive",]
[1] 773
# Eliminar duplicados basados en la columna "sample"
clinical<-clinical[!duplicated(clinical$sample),]
#clinical<-clinical%>%distinct(sample,.keep_all=TRUE)
rownames(clinical)<-clinical$sample
sample_ids<-intersect(rownames(clinical),colnames(mat))
# Filtrar las columnas del conjunto de datos de expresión génica basado en las IDs de las muestras
clinical<-clinical[sample_ids,]
mat<-mat[sample_ids,]
#Normalizacion
mat_norm<-DGEList(mat)
mat_norm<-calcNormFactors(mat_norm,method="TMM")
mat_norm<-cpm(mat_norm,log=TRUE)
mat<-t(mat_norm)
#Anotacion genes
#Creamos la matriz de anotacion con los entrezID y los nombres de gen empleando el paquete AnnotationDbi
entrez<-mapIds(org.Hs.eg.db,keys=colnames(mat),column=c("ENTREZID"),keytipe="SYMBOL")
all(names(entrez)==colnames(mat))
#Eliminamos los no anotados
mat<-mat[,!is.na(entrez)]
entrez<-entrez[!is.na(entrez)]
#La matriz tiene que tener al menos una columna con los entrezId que se llame "EntrezGENE.ID", segun la documentacion de la funcion gene70()
entrezIds<-cbind("EntrezGene.ID"=entrez,"symbol"=colnames(mat))
#Calculo de riesgo con Genefu
#calculamos el riesgo empleando el predictor Mammaprint con la funcion gene 70
risk<-gene70(data=mat,annot=entrezIds,do.mapping=T,verbose=T)
table(risk$risk)
#riesgo binario
#obtenemos los datos de riesgo, siendo 1 riesgo alto y 0 riesgo bajo
clinical$MP_risk<-factor(risk$risk,levels=c(0,1),labels=c("Bajo.riesgo","Alto.riesgo"))
table(clinical$MP_risk)
#Riesgo en modo score
#obtenemos el score de riesgo absoluto
clinical$MP_score<-risk$score
summary(clinical$MP_score)
#Guardamos los datos
write.table(clinical,"clm_clinical_info_TCGA-BRCA_MammaPrintinfo.tsv")
col.names=T,row.names=T,sep=\"t")
data(sig.gene70,package="genefu")
genes_sign<-intersect(colnames(mat),sig.gene70$HUGO.gene.symbol)
mat_filt<-mat[,genes_sign]
saveRDS(mat,file="x.rds")
saveRDS(mat_filt,file="x_filt.rds")
saveRDS(clinical$MP_risk,file="y.rds")
Are you sure you've loaded the library properly? Have you checked if edgeR is ticked in the list of installed packages? Double-check it. If needed, refresh the packages (arrow symbol that rotates)