I have a list of reference genome loci coordinates. The short RNA reads were mapped on these loci by bowtie and samtools.
Instead of simply get the counting number of reads, I hope to generate a matrix like this:
This is just an example showing that read1 aligned two different places in Loci1, 4 different places in loci2, none in loci3.
I had made loci table in a gff format file and loaded it into R. All reads and alignment information is in bam file.
I tried two ways to load bam file into R. One is:
I found the output ignore my reads ID. Only chromosome and positions and other information are available. I need reads ID to count how many different places a read hits within one loci. This seems to be a limitation of the package.
I also tried another method:
This loaded everything and read id was stored in qname. However, I don't know how to perform the further counting to generate the matrix.
Can you give me any suggestion or idea on how to work on it? I tried my best but due to limit knowledge on R and bioconductor, I have suffered from it for several days.
I appreciate if anyone can give me some hints or help me work it out. Thanks!
Alice
Instead of simply get the counting number of reads, I hope to generate a matrix like this:
Code:
read-id Loci1 Loci2 Loci3 ..... 1 2 4 0 2 0 3 4 3 4 2 5 ...
I had made loci table in a gff format file and loaded it into R. All reads and alignment information is in bam file.
I tried two ways to load bam file into R. One is:
Code:
reads<-readBamGappedAlignments("myalign_reads.sorted.bam")
I also tried another method:
Code:
indexBam("myalign_reads.sorted.bam"0 param<-ScanBamParam() reads<-scanBam("myalign_reads.sorted.bam", index="myalign_reads.sorted.bam", param=param)
Can you give me any suggestion or idea on how to work on it? I tried my best but due to limit knowledge on R and bioconductor, I have suffered from it for several days.
I appreciate if anyone can give me some hints or help me work it out. Thanks!
Alice
Comment