How to create automated cbind in a for loop
0
0
Entering edit mode
@arman-shahrisa-7713
Last seen 6.0 years ago

How can I create an automated cbind() to use with for loop?

This is the part that returns error

expression <- cbind(expression, rrM.ProfileData)

The following is my function

genes is a character vector containing names of the desired genes.

cancer is also a character vector string containing names of desired cancers.

csvname will be used as name of exported csv file containing desired expression profile.

makeHeatmapForAllCancers <- function(genes, cancers, csvname, method="mean", zscore.cutoff=2)
{
    print("Make sure the VPN is connected, otherwise the function wont be able to get the required data from http://www.cbioportal.org/")

    ## Type of computing function

    if (method=="mean"){
        cfunc <- function(x) mean(x, na.rm=TRUE)
    } else if (method=="median") {
        cfunc <- function(x) median(x, na.rm=TRUE)
    } else {
        print("method can not be left empety, it should be defined as mean or median")
    }

    ## Getting the required gene expresssion profile

    print("Downloading the required gene expresssion profile from http://www.cbioportal.org/")
    library("cgdsr")
    mycgds = CGDS("http://www.cbioportal.org/")
    colname <- vector("character", length(cancers))
    colname[1] <- "Gene.names"
    samplesize <- vector("numeric", length(cancers))
    cancernames <- vector("character", length(cancers))

    cancers <- (cancers[cancers!="Metastatic Prostate Cancer, SU2C/PCF Dream Team (Robinson et al., Cell 2015)"])
    Metastatic.Prostate.Cancer.mycancerstudy = getCancerStudies(mycgds)[
        which(getCancerStudies(mycgds)[,2]=="Metastatic Prostate Cancer, SU2C/PCF Dream Team (Robinson et al., Cell 2015)"),1]
    Metastatic.Prostate.Cancer.mycaselist = getCaseLists(mycgds,Metastatic.Prostate.Cancer.mycancerstudy)[6,1]
    Metastatic.Prostate.Cancer.mygeneticprofile = getGeneticProfiles(mycgds,Metastatic.Prostate.Cancer.mycancerstudy)[2,1]
    Metastatic.Prostate.Cancer.ProfileData <- t(getProfileData(mycgds,genes,Metastatic.Prostate.Cancer.mygeneticprofile,Metastatic.Prostate.Cancer.mycaselist))
    Expression <- vector("numeric", length=length(genes))
    for(i in 1:length(genes)){
        Expression[i] <- cfunc(as.vector(Metastatic.Prostate.Cancer.ProfileData[i,])[Metastatic.Prostate.Cancer.ProfileData[i,] > zscore.cutoff])
    }
    dim(Expression)  <- c(length(genes),1)

    ## expressionProfile <<- data.frame(expressionProfile, Metastatic.Prostate.Cancer=Expression)

    colname[1] <- "Metastatic.Prostate.Cancer"
    cancernames[1] <- "Metastatic.Prostate.Cancer"
    samplesize[1] <- ncol(Metastatic.Prostate.Cancer.ProfileData)

    ## TCGA

    for(k in 1:length(cancers)){
        mycancerstudy = getCancerStudies(mycgds)[which(getCancerStudies(mycgds)[,2]==as.character(cancers[k])),1]
        mycaselist = getCaseLists(mycgds,mycancerstudy)[which(getCaseLists(mycgds,mycancerstudy)[,2]=="Tumor Samples with mRNA data (RNA Seq V2)"),1]
        mygeneticprofile = getGeneticProfiles(mycgds,mycancerstudy)[which(getGeneticProfiles(mycgds,mycancerstudy)[,2]=="mRNA Expression z-Scores (RNA Seq V2 RSEM)"),1]
        ProfileData <- t(getProfileData(mycgds,genes,mygeneticprofile,mycaselist))
        samplesize[k+2] <- ncol(ProfileData)
        rrM.ProfileData <- vector("numeric", length=length(genes))
        for(j in 1:length(genes)){
            rrM.ProfileData[j] <- cfunc(as.vector(ProfileData[j,])[ProfileData[j,] > zscore.cutoff])
        }

        dim(rrM.ProfileData)  <- c(length(genes),1)
        expression <- cbind(expression, rrM.ProfileData)
        cname <- sapply(strsplit(as.character(cancers[k]), split=" (", fixed=TRUE), function(x) (x[1]))
        colname[k+2] <- as.character(cname)
        cancernames[k+2] <- as.character(cname)
    }

    ## Creating expressionSet class

    print("Creating expressionSet class for obtained data")
    library(Biobase)
    dimnames(expression) <- list(genes, colname)
    expressionpData <- data.frame(Sample.size=samplesize, stringsAsFactors = FALSE)
    colnames(expressionpData) <- colname
    expressionpData <- AnnotatedDataFrame(expressionpData)
    expression.Set <<- ExpressionSet(expression, expressionpData)

    ## Plotting Heatmap

    print("Prepairing Heatmap")
    library(gplots)
    library(RColorBrewer)
    library(rafalib)
    tissue= pData(expression.Set)$tissue
    hmcol <- rev(colorRampPalette(brewer.pal(9, "RdBu"))(100))
    print(heatmap.2(exprs(expression.Set), labCol=tissue, trace="none", symbreaks = T, col=hmcol, cexRow =0.7, cexCol= 0.8, margins= c(17,15)))

    ## Save the expression profile

    genes.list <- t(exprs(expression.Set))
    write.table(genes.list, file=csvname)
    cat("A .csv file entitled", csvname,
        "which contains expression profile for requested genes in all cancers was saved in",
        get.seq(), "directory", sep=" ")

    print("The function was programmed by Arman Shahrisa")
}
loop matix • 2.2k views
ADD COMMENT
1
Entering edit mode

What is the error message you are receiving?

ADD REPLY
0
Entering edit mode

I fixed the problem. It was a "E" instead of "e". I appreciate for your time.

ADD REPLY
0
Entering edit mode

This question has nothing to do with Bioconductor and is of no interest to the Bioconductor community. Can someone simply delete it?

ADD REPLY

Login before adding your answer.

Traffic: 752 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6