1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
|
### =========================================================================
### summarizeVariants methods
### =========================================================================
setMethod("summarizeVariants", c("TxDb", "VCF", "CodingVariants"),
function(query, subject, mode, ...)
{
grl <- cdsBy(query, "tx")
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("TxDb", "VCF", "FiveUTRVariants"),
function(query, subject, mode, ...)
{
grl <- fiveUTRsByTranscript(query)
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("TxDb", "VCF", "ThreeUTRVariants"),
function(query, subject, mode, ...)
{
grl <- threeUTRsByTranscript(query)
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("TxDb", "VCF", "SpliceSiteVariants"),
function(query, subject, mode, ...)
{
grl <- intronsByTranscript(query)
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("TxDb", "VCF", "IntronVariants"),
function(query, subject, mode, ...)
{
grl <- intronsByTranscript(query)
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("TxDb", "VCF", "PromoterVariants"),
function(query, subject, mode, ...)
{
gr <- transcripts(query, columns="tx_id")
grl <- splitAsList(gr, seq_len(length(gr)))
names(grl) <- mcols(gr)$tx_id
callGeneric(grl, subject, mode, ...)
})
setMethod("summarizeVariants", c("GRangesList", "VCF", "VariantType"),
function(query, subject, mode, ...)
{
callGeneric(query, subject, mode=locateVariants, ..., region=mode,
asHits=TRUE)
})
setMethod("summarizeVariants", c("GRangesList", "VCF", "function"),
function(query, subject, mode, ...)
{
if (length(geno(subject)) == 0L) {
warning("No genotypes found in 'query'.")
return(.baseSE(query, subject))
}
## count
hits <- mode(rowRanges(subject), query, ...)
if (length(hits) == 0L)
return(.baseSE(query, subject))
## genotypes
na <- c("0|0", "0/0", "./.", ".|.", ".")
vcf_geno <- geno(subject)$GT[unique(queryHits(hits)), ]
gtype <- as.numeric(!vcf_geno %in% na)
## summarize counts factor-by-sample
fac_x_var <- table(subjectHits(hits), queryHits(hits))
var_x_smp <- matrix(gtype, ncol=ncol(subject))
fac_x_smp <- fac_x_var %*% var_x_smp
rownames(fac_x_smp) <- NULL
SummarizedExperiment(rowRanges=query[unique(subjectHits(hits))],
colData=colData(subject),
metadata=metadata(subject),
assays=SimpleList(counts=fac_x_smp))
})
.baseSE <- function(query, subject, ...)
{
SummarizedExperiment(rowRanges=query, colData=colData(subject),
metadata=metadata(subject),
assays=SimpleList(counts=matrix(NA_integer_,
nrow=length(query), ncol=ncol(subject))))
}
|
