Hi there,

I've just updated Rsamtools to the latest version ("1.22.0") and I've noticed a weird memory leak that is crashing my program.

I have a folder of ~150 bam files that I'm trying to sequentially process, feeding in a data.frame called filter to scanBam(which), that splits the genome into 7744 genomic coordinates.

> head(filter)
        V1       V2       V3                         V4
1 GL896898        0 30180132        GL896898:0-30180132
2 GL896898 30180133 52375790 GL896898:30180133-52375790
3 GL896899        0 24984016        GL896899:0-24984016
4 GL896899 24984017 40953715 GL896899:24984017-40953715
5 GL896905        0 18438707        GL896905:0-18438707
6 GL896905 18438708 27910907 GL896905:18438708-27910907

For some reason, Rsamtools is not releasing RAM after each iteration. I assumed it was my code, but I just ran the following simple test, and unless I'm being completely stupid, Rsamtools is doing something funny. To test this I looped a single bam file multiple time with the following script and monitored RAM usage:

#R implementation of top for system RAM
library(NCmisc)
#Finds total RAM for all objects stored in Rsession
library(pryr)
library(Rsamtools)

filter=read.table('ferret_merged_sce_events.bed')
objectRam=vector()
totalRam=vector()

for(seq in 1:10)
{
     message(paste("iteration", seq))

     #Run scanBam to read in bam file
     bam.dataframe <- scanBam('bamFile.bam', param=ScanBamParam(which=GRanges(seqnames = c(as.character(filter[,1])), ranges = IRanges(c(filter[,2]), c(filter[,3]) )), mapqFilter=10, what=c("rname","pos","strand")))

     #measure RAM across all R objects
     objectRam <- c(objectRam, mem_used()[1]/10^6)
     #measure RAM across system
     totalRam <- c(totalRam,  suppressWarnings(top(CPU=F)$RAM$used*10^6))

     #remove scanBam object to reset for next iteration
     rm(bam.dataframe)

     #empty garbage collection
     gc()

}
plot(objectRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')
plot(totalRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')

This generates the two attached graphs. You'll see that while the RAM used for the R object remains stable, the overall RAM usage of my computer increases, and only decreases when I end my R session. Is there any way of purging the RAM? Is this a bug in Rsamtools or have I done something wrong?

The amount of memory used for all the objects in my R session is constant:

...But the total system usage of RAM increases with each iteration:

I don't think this happened in earlier versions of Rsamtools, as I'd previously analyzed these data without a problem.  I added the rm() and gc() lines in the hope it would fix the problem, but it doesn't seem to change anything. I tried running this without the which= option in scanBamParam() and there doesn't seem to be this memory retention issue, so it's specific for scanBam when which is specified as a parameter. Any help would be greatly appreciated!

EDIT:

I'm adding this to the main section, because it'll probably look choppy in the comments, but I changed things around as recommended by Martin Morgan, and for some reason, I'm only getting an increase in system memory if param is specified within the loop, but not if it's specified outside of it.

What I was doing was the following, which is eating my system memory:

library(NCmisc)
library(pryr)
library(Rsamtools)
library(RNAseqData.HNRNPC.bam.chr14)
options(scipen=20)
#make equivalent table of 10,000 regions in data
filter=data.frame(chr='chr14', start=seq(0, 105000000, 10000), end=seq(10000, 105010000, 10000))
objectRam=vector()
totalRam=vector()
for(fileName in RNAseqData.HNRNPC.bam.chr14_BAMFILES)
{
    bam.dataframe <- scanBam(fileName, param=ScanBamParam(which=GRanges(seqnames=c(as.character(filter[,1])), ranges=IRanges(c(filter[,2]), c(filter[,3]) )), mapqFilter=10, what=c("rname","pos","strand")))
    objectRam <- c(objectRam, mem_used()[1]/10^6)
    totalRam <- c(totalRam,  suppressWarnings(top(CPU=F)$RAM$used*10^6))
    rm(bam.dataframe)
    gc()
}
plot(objectRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')
plot(totalRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')

but if I just take param out of the loop, my memory issues seem to go away:

library(NCmisc)
library(pryr)
library(Rsamtools)
library(RNAseqData.HNRNPC.bam.chr14)

options(scipen=20)
#make equivalent table of 10,000 regions in data
filter=data.frame(chr='chr14', start=seq(0, 105000000, 10000), end=seq(10000, 105010000, 10000))
objectRam=vector()
totalRam=vector()

parameters=ScanBamParam(which=GRanges(seqnames=c(as.character(filter[,1])), ranges=IRanges(c(filter[,2]), c(filter[,3]) )), mapqFilter=10, what=c("rname","pos","strand"))

for(fileName in RNAseqData.HNRNPC.bam.chr14_BAMFILES)
{
    bam.dataframe <- scanBam(fileName, param=parameters)
    objectRam <- c(objectRam, mem_used()[1]/10^6)
    totalRam <- c(totalRam,  suppressWarnings(top(CPU=F)$RAM$used*10^6))
    rm(bam.dataframe)
    gc()
}

plot(objectRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')
plot(totalRam, type='l', xlab='iteration', ylab='Ram usage (Mb)')

I'm not sure why this is, but it seems to be working much better, so I'll try it on my actual data.  Thanks so much!

PS: I didn't know message() natively read stored variables. Thanks! That'll clean things up!