translationalEff {ORFik} | R Documentation |
Uses RnaSeq and RiboSeq to get translational efficiency of every element in 'grl'. Translational efficiency is defined as:
(density of RPF within ORF) / (RNA expression of ORFs transcript)
translationalEff(grl, RNA, RFP, tx, with.fpkm = FALSE, pseudoCount = 0)
grl |
a |
RNA |
RnaSeq reads as GAlignment, GRanges or GRangesList object |
RFP |
RiboSeq reads as GAlignment, GRanges or GRangesList object |
tx |
a GRangesList of the transcripts. If you used cage data, then the tss for the the leaders have changed, therefor the tx lengths have changed. To account for that call: ' translationalEff(grl, RNA, RFP, tx = extendLeaders(tx, cageFiveUTRs)) ' where cageFiveUTRs are the reannotated by CageSeq data leaders. |
with.fpkm |
logical F, if true return the fpkm values together with translational efficiency |
pseudoCount |
an integer, 0, set it to 1 if you want to avoid NA and inf values. It also helps against bias from low depth libraries. |
a numeric vector of fpkm ratios, if with.fpkm is TRUE, return a data.table with te and fpkm values
doi: 10.1126/science.1168978
Other features: computeFeaturesCage
,
computeFeatures
,
disengagementScore
,
distToCds
, distToTSS
,
entropy
, floss
,
fpkm_calc
, fpkm
,
fractionLength
,
initiationScore
,
insideOutsideORF
, isInFrame
,
isOverlapping
,
kozakSequenceScore
, orfScore
,
rankOrder
,
ribosomeReleaseScore
,
ribosomeStallingScore
,
startRegionCoverage
,
startRegion
, subsetCoverage
ORF <- GRanges(seqnames = "1", ranges = IRanges(start = c(1, 10, 20), end = c(5, 15, 25)), strand = "+") grl <- GRangesList(tx1_1 = ORF) RFP <- GRanges("1", IRanges(25, 25), "+") RNA <- GRanges("1", IRanges(1, 50), "+") tx <- GRangesList(tx1 = GRanges("1", IRanges(1, 50), "+")) # grl must have same names as cds + _1 etc, so that they can be matched. te <- translationalEff(grl, RNA, RFP, tx, with.fpkm = TRUE, pseudoCount = 1) te$fpkmRFP te$te