Fit the standard mass action kinetics model of gene expression by different methods. Some methods require steady state assumptions, for others data must be properly normalized. The parameters are fit per Condition.
Usage
FitKinetics(
data,
name.prefix = "kinetics",
type = c("nlls", "ntr", "lm", "chase"),
slot = DefaultSlot(data),
time = Design$dur.4sU,
CI.size = 0.95,
return.fields = c("Synthesis", "Half-life"),
return.extra = NULL,
...
)
Arguments
- data
A grandR object
- name.prefix
the prefix of the analysis name to be stored in the grandR object
- type
Which method to use (either one of "full","ntr","lm", "chase")
- slot
The data slot to take expression values from
- time
The column in the column annotation table representing the labeling duration
- CI.size
A number between 0 and 1 representing the size of the confidence interval
- return.fields
which statistics to return (see details)
- return.extra
additional statistics to return (see details)
- ...
forwarded to
FitKineticsGeneNtr
,FitKineticsGeneLeastSquares
orFitKineticsGeneLogSpaceLinear
Details
The start of labeling for all samples should be the same experimental time point. The fit gets more precise with multiple samples from multiple labeling durations.
The standard mass action kinetics model of gene expression arises from the following differential equation:
$$df/dt = s - d f(t)$$
This model assumes constant synthesis and degradation rates. Based on this, there are different ways for fitting the parameters:
FitKineticsGeneLeastSquares: non-linear least squares fit on the full model; depends on proper normalization; can work without steady state; assumption of homoscedastic gaussian errors is theoretically not justified
FitKineticsGeneLogSpaceLinear: linear model fit on the old RNA; depends on proper normalization; assumes steady state for estimating the synthesis rate; assumption of homoscedastic gaussian errors in log space is problematic and theoretically not justified
FitKineticsGeneNtr: maximum a posteriori fit on the NTR posterior transformed to the degradation rate; as it is based on the NTR only, it is independent on proper normalization; assumes steady state; theoretically well justified
Pulse-chase designs are fit using FitKineticsGeneLeastSquares while only considering the drop of labeled RNA. Note that in this case the notion "new" / "old" RNA is misleading, since labeled RNA corresponds to pre-existing RNA!
This function is flexible in what to put in the analysis table. You can specify the statistics using return.fields and return.extra (see kinetics2vector
)
Examples
sars <- ReadGRAND(system.file("extdata", "sars.tsv.gz", package = "grandR"),
design=c("Cell",Design$dur.4sU,Design$Replicate))
#> Warning: Duplicate gene symbols (n=1, e.g. MATR3) present, making unique!
sars <- FilterGenes(sars,use=1:10)
sars<-FitKinetics(sars,name="kinetics.ntr",type='ntr')
sars<-Normalize(sars)
sars<-FitKinetics(sars,name="kinetics.nlls",type='nlls')
sars<-FitKinetics(sars,name="kinetics.lm",type='lm')
head(GetAnalysisTable(sars,columns="Half-life"))
#> Gene Symbol Length Type kinetics.ntr.Half-life
#> UHMK1 ENSG00000152332 UHMK1 8478 Cellular 7.2504385
#> ATF3 ENSG00000162772 ATF3 2103 Cellular 0.7928108
#> PABPC4 ENSG00000090621 PABPC4 3592 Cellular 5.3823958
#> ROR1 ENSG00000185483 ROR1 5832 Cellular 3.4115435
#> ZC3H11A ENSG00000058673 ZC3H11A 11825 Cellular 2.2714713
#> ZBED6 ENSG00000257315 ZBED6 12481 Cellular 2.2385003
#> kinetics.nlls.Half-life kinetics.lm.Half-life
#> UHMK1 4.822883 5.5196842
#> ATF3 0.860518 0.8955509
#> PABPC4 3.792914 5.4373019
#> ROR1 2.278106 1.1723454
#> ZC3H11A 1.728668 0.9935412
#> ZBED6 1.694603 0.9029407