Skip to content

Commit

Permalink
Add unified pair cache across bam/cram
Browse files Browse the repository at this point in the history
  • Loading branch information
@cmdcolin
cmdcolin committed Nov 18, 2018
1 parent c20d970 commit 1d893b2
Show file tree
Hide file tree
Showing 3 changed files with 16 additions and 84 deletions.
21 changes: 1 addition & 20 deletions src/JBrowse/Store/SeqFeature/BAM.js
View file
Expand Up @@ -150,7 +150,6 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, In
})

this.storeTimeout = args.storeTimeout || 3000;
this.featureCache = {}
},

// process the parsed SAM header from the bam file
Expand Down Expand Up @@ -220,7 +219,7 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, In
.then(records => {
if(query.viewAsPairs) {
const recs = records.map(f => this._bamRecordToFeature(f))
this.pairReads(records, featCallback, endCallback, errorCallback)
this.pairFeatures(query, recs, featCallback, endCallback, errorCallback)
} else {
for(let i = 0; i < records.length; i++) {
featCallback(this._bamRecordToFeature(records[i]))
Expand All @@ -230,24 +229,6 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, In
}).catch(errorCallback)
},

cleanFeatureCache(query) {
Object.entries(this.featureCache).forEach(([k, v]) => {
if((v._get('end') < query.start) || (v._get('start') > query.end)) {
delete this.featureCache[k]
}
})
},

getStatsForPairCache() {
if(Object.keys(this.featureCache).length > 400) {
var tlens = Object.entries(this.featureCache).map(([k, v]) => Math.abs(v.get('template_length'))).filter(x => x < MAX_INSERT_SIZE_FOR_STATS).sort((a, b) => a - b)
return {
upper: Util.percentile(tlens, 0.995),
lower: Util.percentile(tlens, 0.005)
}
}
return { upper: Infinity, lower: 0 }
},

_bamRecordToFeature(record) {
return new BamSlightlyLazyFeature(record, this)
Expand Down
57 changes: 4 additions & 53 deletions src/JBrowse/Store/SeqFeature/CRAM.js
View file
@@ -1,36 +1,11 @@
const LRU = cjsRequire('lru-cache')
const { IndexedCramFile, CraiIndex } = cjsRequire('@gmod/cram')

const { Buffer } = cjsRequire('buffer')

const cramIndexedFilesCache = LRU(5)

const BlobFilehandleWrapper = cjsRequire('../../Model/BlobFilehandleWrapper')
const MAX_INSERT_SIZE_FOR_STATS = 30000

class PairedCramRead {
id() {
return Math.min(this.read1.id(), this.read2.id())
}
get(field) {
return this._get(field.toLowerCase())
}
_get(field) {
if(field === 'start') {
return Math.min(this.read1._get('start'), this.read2._get('start'))
} else if(field === 'end') {
return Math.max(this.read1._get('end'), this.read2._get('end'))
} else if(field === 'name') {
return this.read1._get('name')
} else if(field === 'pair_orientation') {
return this.read1._get('pair_orientation')
} else if(field === 'template_length') {
return this.read1._get('template_length')
}
}
pairedFeature() { return true }
children() {}
}
class CramSlightlyLazyFeature {

_get_name() { return this.record.readName }
Expand Down Expand Up @@ -100,13 +75,6 @@ class CramSlightlyLazyFeature {
pairedFeature() { return false }
}

function canBePaired(alignment) {
return alignment.isPaired() &&
!alignment.isMateUnmapped() &&
alignment.sequenceId === alignment.mate.sequenceId &&
(alignment.isRead1() || alignment.isRead2()) &&
!(alignment.isSecondary() || alignment.isSupplementary());
}

define( [
'dojo/_base/declare',
Expand All @@ -116,6 +84,7 @@ define( [
'JBrowse/Store/DeferredStatsMixin',
'JBrowse/Store/DeferredFeaturesMixin',
'JBrowse/Store/SeqFeature/GlobalStatsEstimationMixin',
'JBrowse/Store/SeqFeature/_PairCache',
'JBrowse/Model/XHRBlob',
'JBrowse/Model/SimpleFeature',
],
Expand All @@ -127,11 +96,12 @@ define( [
DeferredStatsMixin,
DeferredFeaturesMixin,
GlobalStatsEstimationMixin,
PairCache,
XHRBlob,
SimpleFeature,
) {

return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, GlobalStatsEstimationMixin ],
return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, GlobalStatsEstimationMixin, PairCache ],

/**
* @lends JBrowse.Store.SeqFeature.CRAM
Expand Down Expand Up @@ -199,7 +169,6 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, Gl
})

this.storeTimeout = args.storeTimeout || 3000;
this.featureCache = {}
},

// process the parsed SAM header from the cram file
Expand Down Expand Up @@ -325,7 +294,7 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, Gl
.then(records => {
if(query.viewAsPairs) {
const recs = records.map(f => this._cramRecordToFeature(f))
this.pairReads(records, featCallback, endCallback, errorCallback)
this.pairFeatures(query, recs, featCallback, endCallback, errorCallback)
} else {
for(let i = 0; i < records.length; i++) {
featCallback(this._cramRecordToFeature(records[i]))
Expand All @@ -341,24 +310,6 @@ return declare( [ SeqFeatureStore, DeferredStatsMixin, DeferredFeaturesMixin, Gl
})
},

cleanFeatureCache(query) {
Object.entries(this.featureCache).forEach(([k, v]) => {
if((v.get('end') < query.start) || (v.get('start') > query.end)) {
delete this.featureCache[k]
}
})
},
getStatsForPairCache() {
if(Object.keys(this.featureCache).length > 400) {
var tlens = Object.entries(this.featureCache).map(([k, v]) => Math.abs(v.get('template_length'))).filter(x => x < MAX_INSERT_SIZE_FOR_STATS).sort((a, b) => a - b)
return {
upper: Util.percentile(tlens, 0.995),
lower: Util.percentile(tlens, 0.005)
}
}
return { upper: Infinity, lower: 0 }
},

_cramRecordToFeature(record) {
return new CramSlightlyLazyFeature(record, this)
},
Expand Down
22 changes: 11 additions & 11 deletions src/JBrowse/Store/SeqFeature/_PairCache.js
View file
Expand Up @@ -23,11 +23,11 @@ class PairedRead {
}

function canBePaired(alignment) {
return alignment.isPaired() &&
!alignment.isMateUnmapped() &&
alignment._refID === alignment._next_refid() &&
(alignment.isRead1() || alignment.isRead2()) &&
!(alignment.isSecondary() || alignment.isSupplementary());
return alignment.get('multi_segment_template') &&
!alignment.get('multi_segment_next_segment_unmapped') &&
alignment.get('seq_id') === alignment.get('next_seq_id') &&
(alignment.get('multi_segment_first') || alignment.get('multi_segment_last')) &&
!(alignment.get('secondary_alignment') || alignment.get('supplementary_alignment'))
}

define( [
Expand All @@ -52,17 +52,17 @@ return declare(null, {


// called by getFeatures from the DeferredFeaturesMixin
pairFeatures(features, featCallback, endCallback, errorCallback, featTransform) {
pairFeatures(query, records, featCallback, endCallback, errorCallback, featTransform) {
const pairCache = {};
for(let i = 0; i < records.length; i++) {
let feat
if (canBePaired(records[i])) {
let name = records[i]._get('name')
feat = pairCache[name]
if (feat) {
if(records[i].isRead1()) {
if(records[i].get('multi_segment_first')) {
feat.read1 = records[i]
} else if(records[i].isRead2()) {
} else if(records[i].get('multi_segment_last')) {
feat.read2 = records[i]
} else {
console.log('unable to pair read',records[i])
Expand All @@ -73,10 +73,10 @@ return declare(null, {
}
}
else {
feat = new PairedBamRead()
if(records[i].isRead1()) {
feat = new PairedRead()
if(records[i].get('multi_segment_first')) {
feat.read1 = records[i]
} else if(records[i].isRead2()) {
} else if(records[i].get('multi_segment_last')) {
feat.read2 = records[i]
} else {
console.log('unable to pair read', records[i])
Expand Down

0 comments on commit 1d893b2

Please sign in to comment.