///|
pub fn Annotation::build_interval_index(self : Annotation) -> IntervalIndex {
{ features: self.features }
}
///|
pub fn IntervalIndex::query(
self : IntervalIndex,
seqid~ : String,
start~ : Int,
end~ : Int,
) -> Array[IntervalHit] {
let hits : Array[IntervalHit] = []
for index, feature in self.features {
if feature.overlaps(seqid~, start~, end~) {
hits.push({ feature, index })
}
}
hits
}
///|
pub fn IntervalIndex::query_span(
self : IntervalIndex,
span : Span,
) -> Array[IntervalHit] {
self.query(seqid=span.seqid, start=span.start, end=span.end)
}
///|
pub fn IntervalIndex::count_seqid(self : IntervalIndex, seqid : String) -> Int {
let mut count = 0
for feature in self.features {
if feature.seqid == seqid {
count += 1
}
}
count
}
///|
pub fn IntervalIndex::seqids(self : IntervalIndex) -> Array[String] {
let result : Array[String] = []
for feature in self.features {
if !result.contains(feature.seqid) {
result.push(feature.seqid)
}
}
result
}
///|
pub fn IntervalIndex::features_on(
self : IntervalIndex,
seqid : String,
) -> Array[Feature] {
let result : Array[Feature] = []
for feature in self.features {
if feature.seqid == seqid {
result.push(feature)
}
}
result
}
///|
fn transcript_name(tx : TranscriptModel) -> String {
tx.name.unwrap_or(tx.id)
}
///|
fn bed12_blocks(tx : TranscriptModel) -> (Array[String], Array[String]) {
let sizes : Array[String] = []
let starts : Array[String] = []
let tx_start = tx.feature.start
for exon in tx.exons {
sizes.push("\{exon.length()}")
starts.push("\{exon.start - tx_start}")
}
(sizes, starts)
}
///|
fn transcript_cds_bounds(tx : TranscriptModel) -> (Int, Int) {
if tx.cds.length() == 0 {
(tx.feature.to_zero_based_start(), tx.feature.to_half_open_end())
} else {
let mut start = tx.cds[0].start
let mut end = tx.cds[0].end
for cds in tx.cds {
if cds.start < start {
start = cds.start
}
if cds.end > end {
end = cds.end
}
}
(start - 1, end)
}
}
///|
pub fn TranscriptModel::to_bed12(self : TranscriptModel) -> Bed12Record {
let (sizes, starts) = bed12_blocks(self)
let (thick_start, thick_end) = transcript_cds_bounds(self)
{
chrom: self.feature.seqid,
chrom_start: self.feature.to_zero_based_start(),
chrom_end: self.feature.to_half_open_end(),
name: transcript_name(self),
score: "0",
strand: self.feature.strand.to_string(),
thick_start,
thick_end,
item_rgb: "0",
block_count: sizes.length(),
block_sizes: sizes,
block_starts: starts,
}
}
///|
pub fn Bed12Record::to_line(self : Bed12Record) -> String {
[
self.chrom,
"\{self.chrom_start}",
"\{self.chrom_end}",
self.name,
self.score,
self.strand,
"\{self.thick_start}",
"\{self.thick_end}",
self.item_rgb,
"\{self.block_count}",
self.block_sizes.join(","),
self.block_starts.join(","),
].join("\t")
}
///|
pub fn Annotation::to_bed12(self : Annotation) -> String {
let lines : Array[String] = []
for gene in self.to_gene_models() {
for tx in gene.transcripts {
lines.push(tx.to_bed12().to_line())
}
}
if lines.length() == 0 {
""
} else {
lines.join("\n") + "\n"
}
}