Lack-of-fit test for models fitted to data with replicates

+ + +

This is a generic function with a method currently only defined for mkinfit +objects. It fits an anova model to the data contained in the object and +compares the likelihoods using the likelihood ratio test +lrtest.default from the lmtest package.

+ +

loftest(object, ...)
+
+# S3 method for mkinfit
+loftest(object, ...)

+ +

Arguments

+ + + + + + + + + + +

object	A model object with a defined loftest method
...	Not used

+ +

Details

+ +

The anova model is interpreted as the simplest form of an mkinfit model, +assuming only a constant variance about the means, but not enforcing any +structure of the means, so we have one model parameter for every mean +of replicate samples.

Examples

# \dontrun{
+test_data <- subset(synthetic_data_for_UBA_2014[[12]]$data, name == "parent")
+sfo_fit <- mkinfit("SFO", test_data, quiet = TRUE)
+plot_res(sfo_fit) # We see a clear pattern in the residuals
loftest(sfo_fit)  # We have a clear lack of fit
#> Likelihood ratio test
+#> 
+#> Model 1: ANOVA with error model const
+#> Model 2: SFO with error model const
+#>   #Df  LogLik Df  Chisq Pr(>Chisq)    
+#> 1  10 -40.710                         
+#> 2   3 -63.954 -7 46.487  7.027e-08 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
+# We try a different model (the one that was used to generate the data)
+dfop_fit <- mkinfit("DFOP", test_data, quiet = TRUE)
+plot_res(dfop_fit) # We don't see systematic deviations, but heteroscedastic residuals
# therefore we should consider adapting the error model, although we have
+loftest(dfop_fit) # no lack of fit
#> Likelihood ratio test
+#> 
+#> Model 1: ANOVA with error model const
+#> Model 2: DFOP with error model const
+#>   #Df  LogLik Df Chisq Pr(>Chisq)
+#> 1  10 -40.710                    
+#> 2   5 -42.453 -5 3.485     0.6257
#
+# This is the anova model used internally for the comparison
+test_data_anova <- test_data
+test_data_anova$time <- as.factor(test_data_anova$time)
+anova_fit <- lm(value ~ time, data = test_data_anova)
+summary(anova_fit)
#> 
+#> Call:
+#> lm(formula = value ~ time, data = test_data_anova)
+#> 
+#> Residuals:
+#>     Min      1Q  Median      3Q     Max 
+#> -6.1000 -0.5625  0.0000  0.5625  6.1000 
+#> 
+#> Coefficients:
+#>             Estimate Std. Error t value Pr(>|t|)    
+#> (Intercept)  103.150      2.323  44.409 7.44e-12 ***
+#> time1        -19.950      3.285  -6.073 0.000185 ***
+#> time3        -50.800      3.285 -15.465 8.65e-08 ***
+#> time7        -68.500      3.285 -20.854 6.28e-09 ***
+#> time14       -79.750      3.285 -24.278 1.63e-09 ***
+#> time28       -86.000      3.285 -26.181 8.35e-10 ***
+#> time60       -94.900      3.285 -28.891 3.48e-10 ***
+#> time90       -98.500      3.285 -29.986 2.49e-10 ***
+#> time120     -100.450      3.285 -30.580 2.09e-10 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+#> 
+#> Residual standard error: 3.285 on 9 degrees of freedom
+#> Multiple R-squared:  0.9953,	Adjusted R-squared:  0.9912 
+#> F-statistic: 240.5 on 8 and 9 DF,  p-value: 1.417e-09
+#> 
logLik(anova_fit) # We get the same likelihood and degrees of freedom
#> 'log Lik.' -40.71015 (df=10)
#
+test_data_2 <- synthetic_data_for_UBA_2014[[12]]$data
+m_synth_SFO_lin <- mkinmod(parent = list(type = "SFO", to = "M1"),
+  M1 = list(type = "SFO", to = "M2"),
+  M2 = list(type = "SFO"), use_of_ff = "max")
#> Successfully compiled differential equation model from auto-generated C code.
sfo_lin_fit <- mkinfit(m_synth_SFO_lin, test_data_2, quiet = TRUE)
+plot_res(sfo_lin_fit) # not a good model, we try parallel formation
loftest(sfo_lin_fit)
#> Likelihood ratio test
+#> 
+#> Model 1: ANOVA with error model const
+#> Model 2: m_synth_SFO_lin with error model const and fixed parameter(s) M1_0, M2_0
+#>   #Df   LogLik  Df  Chisq Pr(>Chisq)    
+#> 1  28  -93.606                          
+#> 2   7 -171.927 -21 156.64  < 2.2e-16 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
+m_synth_SFO_par <- mkinmod(parent = list(type = "SFO", to = c("M1", "M2")),
+  M1 = list(type = "SFO"),
+  M2 = list(type = "SFO"), use_of_ff = "max")
#> Successfully compiled differential equation model from auto-generated C code.
sfo_par_fit <- mkinfit(m_synth_SFO_par, test_data_2, quiet = TRUE)
+plot_res(sfo_par_fit) # much better for metabolites
loftest(sfo_par_fit)
#> Likelihood ratio test
+#> 
+#> Model 1: ANOVA with error model const
+#> Model 2: m_synth_SFO_par with error model const and fixed parameter(s) M1_0, M2_0
+#>   #Df   LogLik  Df  Chisq Pr(>Chisq)    
+#> 1  28  -93.606                          
+#> 2   7 -156.331 -21 125.45  < 2.2e-16 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
+m_synth_DFOP_par <- mkinmod(parent = list(type = "DFOP", to = c("M1", "M2")),
+  M1 = list(type = "SFO"),
+  M2 = list(type = "SFO"), use_of_ff = "max")
#> Successfully compiled differential equation model from auto-generated C code.
dfop_par_fit <- mkinfit(m_synth_DFOP_par, test_data_2, quiet = TRUE)
+plot_res(dfop_par_fit) # No visual lack of fit
loftest(dfop_par_fit)  # no lack of fit found by the test
#> Likelihood ratio test
+#> 
+#> Model 1: ANOVA with error model const
+#> Model 2: m_synth_DFOP_par with error model const and fixed parameter(s) M1_0, M2_0
+#>   #Df   LogLik  Df  Chisq Pr(>Chisq)
+#> 1  28  -93.606                      
+#> 2   9 -102.763 -19 18.313     0.5016
#
+# The anova model used for comparison in the case of transformation products
+test_data_anova_2 <- dfop_par_fit$data
+test_data_anova_2$variable <- as.factor(test_data_anova_2$variable)
+test_data_anova_2$time <- as.factor(test_data_anova_2$time)
+anova_fit_2 <- lm(observed ~ time:variable - 1, data = test_data_anova_2)
+summary(anova_fit_2)
#> 
+#> Call:
+#> lm(formula = observed ~ time:variable - 1, data = test_data_anova_2)
+#> 
+#> Residuals:
+#>     Min      1Q  Median      3Q     Max 
+#> -6.1000 -0.5875  0.0000  0.5875  6.1000 
+#> 
+#> Coefficients: (2 not defined because of singularities)
+#>                        Estimate Std. Error t value Pr(>|t|)    
+#> time0:variableparent    103.150      1.573  65.562  < 2e-16 ***
+#> time1:variableparent     83.200      1.573  52.882  < 2e-16 ***
+#> time3:variableparent     52.350      1.573  33.274  < 2e-16 ***
+#> time7:variableparent     34.650      1.573  22.024  < 2e-16 ***
+#> time14:variableparent    23.400      1.573  14.873 6.35e-14 ***
+#> time28:variableparent    17.150      1.573  10.901 5.47e-11 ***
+#> time60:variableparent     8.250      1.573   5.244 1.99e-05 ***
+#> time90:variableparent     4.650      1.573   2.956 0.006717 ** 
+#> time120:variableparent    2.700      1.573   1.716 0.098507 .  
+#> time0:variableM1             NA         NA      NA       NA    
+#> time1:variableM1         11.850      1.573   7.532 6.93e-08 ***
+#> time3:variableM1         22.700      1.573  14.428 1.26e-13 ***
+#> time7:variableM1         33.050      1.573  21.007  < 2e-16 ***
+#> time14:variableM1        31.250      1.573  19.863  < 2e-16 ***
+#> time28:variableM1        18.900      1.573  12.013 7.02e-12 ***
+#> time60:variableM1         7.550      1.573   4.799 6.28e-05 ***
+#> time90:variableM1         3.850      1.573   2.447 0.021772 *  
+#> time120:variableM1        2.050      1.573   1.303 0.204454    
+#> time0:variableM2             NA         NA      NA       NA    
+#> time1:variableM2          6.700      1.573   4.259 0.000254 ***
+#> time3:variableM2         16.750      1.573  10.646 8.93e-11 ***
+#> time7:variableM2         25.800      1.573  16.399 6.89e-15 ***
+#> time14:variableM2        28.600      1.573  18.178 6.35e-16 ***
+#> time28:variableM2        25.400      1.573  16.144 9.85e-15 ***
+#> time60:variableM2        21.600      1.573  13.729 3.81e-13 ***
+#> time90:variableM2        17.800      1.573  11.314 2.51e-11 ***
+#> time120:variableM2       14.100      1.573   8.962 2.79e-09 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+#> 
+#> Residual standard error: 2.225 on 25 degrees of freedom
+#> Multiple R-squared:  0.9979,	Adjusted R-squared:  0.9957 
+#> F-statistic: 469.2 on 25 and 25 DF,  p-value: < 2.2e-16
+#> 
# }
+

+ +

object	An `mkinfit` object	An `mkinfit` object, or an `mmkin` column +object containing two fits to the same data.
object_2

object	An mmkin column object, containing two or more +`mkinfit` models that have been fitted to the same data, +or an mkinfit object. In the latter case, further mkinfit +objects fitted to the same data should be specified +as dots arguments.
...	Not used in the method for mmkin column objects, +further mkinfit objects in the method for mkinfit objects.

-- cgit v1.2.1 From 49d80b149ec36eb967cb6da3a8b6726e48bf0ecc Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Tue, 10 Dec 2019 20:46:05 +0100 Subject: Update and amend the README --- docs/index.html | 22 ++++++++++++++++++++-- 1 file changed, 20 insertions(+), 2 deletions(-) (limited to 'docs') diff --git a/docs/index.html b/docs/index.html index c89fc122..2e960df1 100644 --- a/docs/index.html +++ b/docs/index.html @@ -151,8 +151,8 @@

Credits and historical remarks

mkin would not be possible without the underlying software stack consisting of R and the packages deSolve and FME, to say the least.

It could not have been written without me being introduced to regulatory fate modelling of pesticides by Adrian Gurney during my time at Harlan Laboratories Ltd (formerly RCC Ltd). mkin greatly profits from and largely follows the work done by the FOCUS Degradation Kinetics Workgroup, as detailed in their guidance document from 2006, slightly updated in 2011 and in 2014.

mkin would not be possible without the underlying software stack consisting of R and the package deSolve. In previous version, mkin was also using the functionality of the FME package.

mkin could not have been written without me being introduced to regulatory fate modelling of pesticides by Adrian Gurney during my time at Harlan Laboratories Ltd (formerly RCC Ltd). mkin greatly profits from and largely follows the work done by the FOCUS Degradation Kinetics Workgroup, as detailed in their guidance document from 2006, slightly updated in 2011 and in 2014.

Also, it was inspired by the first version of KinGUI developed by BayerCropScience, which is based on the MatLab runtime environment.

The companion package kinfit (now deprecated) was started in 2008 and first published on CRAN on 01 May 2010.

The first mkin code was published on 11 May 2010 and the first CRAN version on 18 May 2010.

@@ -160,6 +160,24 @@

Somewhat in parallel, Syngenta has sponsored the development of an mkin and KinGUII based GUI application called CAKE, which also adds IRLS and MCMC, is more limited in the model formulation, but puts more weight on usability. CAKE is available for download from the CAKE website, where you can also find a zip archive of the R scripts derived from mkin, published under the GPL license.

Finally, there is KineticEval, which contains a further development of the scripts used for KinGUII, so the different tools will hopefully be able to learn from each other in the future as well.

+References

+ + + +

Ranke J, Meinecke S (2019) + Error Models for the Kinetic Evaluation of Chemical Degradation Data + Environments + 6 (12) 124 + doi:10.3390/environments6120124 +

Ranke J, Wöltjen J, Meinecke S (2018) + Comparison of software tools for kinetic evaluation of chemical degradation data + Environmental Sciences Europe + 30 17 + doi:10.1186/s12302-018-0145-1 +

Development

-- cgit v1.2.1 From 7c8e06d4bae9ec80574147005b7e933937447220 Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Tue, 10 Dec 2019 20:56:10 +0100 Subject: More additions to the README --- docs/index.html | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'docs') diff --git a/docs/index.html b/docs/index.html index 2e960df1..fcb9c90b 100644 --- a/docs/index.html +++ b/docs/index.html @@ -151,7 +151,7 @@

Credits and historical remarks

mkin would not be possible without the underlying software stack consisting of R and the package deSolve. In previous version, mkin was also using the functionality of the FME package.

mkin would not be possible without the underlying software stack consisting of, among others, R and the package deSolve. In previous version, mkin was also using the functionality of the FME package. Please refer to the package page on CRAN for the full list of imported and suggested R packages. Also, Debian Linux, the vim editor and the Nvim-R plugin have been invaluable in its development.

Also, it was inspired by the first version of KinGUI developed by BayerCropScience, which is based on the MatLab runtime environment.

The companion package kinfit (now deprecated) was started in 2008 and first published on CRAN on 01 May 2010.

-- cgit v1.2.1 From 1868c1c6b98afa4c8a11b7c065d717bfb4ec1a8e Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Mon, 16 Dec 2019 02:51:19 +0100 Subject: Argument ymax for plot.mmkin --- docs/reference/plot.mmkin.html | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) (limited to 'docs') diff --git a/docs/reference/plot.mmkin.html b/docs/reference/plot.mmkin.html index 8b68cfae..18907aa2 100644 --- a/docs/reference/plot.mmkin.html +++ b/docs/reference/plot.mmkin.html @@ -73,7 +73,7 @@ the fit of at least one model to the same dataset is shown." />

@@ -147,7 +147,7 @@ the fit of at least one model to the same dataset is shown.

plot(x, main = "auto", legends = 1, resplot = c("time", "errmod"), show_errmin = TRUE, errmin_var = "All data", errmin_digits = 3, cex = 0.7, - rel.height.middle = 0.9, ...)

Lack-of-fit test for models fitted to data with replicates

Arguments

Details

See also

Examples

Contents

Arguments

Calculate Akaike weights for model averaging

Arguments

References

Examples

Contents

Development

References

Examples

Two component error model

Two-component error model

Credits and historical remarks

+References

Development

Credits and historical remarks

Arguments