* Apologies for multiple posting *
I attached to my previous e-mail a .r file, and it was not permitted by the
rules of the mailing lis. Again, please receive my sincere apologies for
this.
I re-send again the e-mail with .txt attachemnt in the hope someone an help
me to solve my problem.
I'm trying to fit a set an ODE to an experimental time series. In the
attachment you find the R code I wrote using modFit and modCost of FME
package and the file of the time series.
When I run summary(Fit) I obtain this error message, and the values of the
parameters are equal to the initial guesses I gave to them.
The problem is not due to the fact that I have only one equation (I tried
also with more equations, but I still obtain this error).
I would appreciate if someone could help me in understanding the reason of
the error and in fixing it.
Thanks for your attention,
Paola Lecca.
Here the error:
> summary(Fit)
Parameters:
Estimate Std. Error t value Pr(>|t|)
pro1_strength 1 NA NA NA
Residual standard error: 2.124 on 10 degrees of freedom
Error in cov2cor(x$cov.unscaled) : 'V' is not a square numeric matrix
In addition: Warning message:
In summary.modFit(Fit) : Cannot estimate covariance; system is singular
______________________________________________
R-help at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
*Paola Lecca, PhD*
*The Microsoft Research - University of Trento*
*Centre for Computational and Systems Biology*
*Piazza Manci 17 38123 Povo/Trento, Italy*
*Phome: +39 0461282843*
*Fax: +39 0461282814*
-------------- next part --------------
time pp1_mrna
0 0
2 2.754
4 2.958
6 4.058
8 3.41
10 3.459
12 2.453
14 1.234
16 2.385
18 3.691
20 3.252
-------------- next part --------------
require(deSolve)
require(FME)
##################################################################################
# PART 1
#
##################################################################################
# Differential equations
model_1_part_1 <- function(t, S, parameters)
{
with(as.list(parameters), {
#
cod1 = pro1_strength
#
pp1_mrna_degradation_rate <- 1
###############################################################
#
v1 = cod1
v2 = pp1_mrna_degradation_rate * S[1]
#
#################################################################
#
dS1 = v1 - v2
#
####################################################################
#
list(c(dS1))
})
}
# Parameters
parms_part_1 <- c(pro1_strength = 1.0)
# Initial values of the species concentration
S <- c(pp1_mrna = 0)
times <- seq(0, 20, by = 2)
# Solve the system
ode_solutions_part_1 <- ode(S, times, model_1_part_1, parms = parms_part_1)
ode_solutions_part_1
summary(ode_solutions_part_1)
## Default plot method
plot(ode_solutions_part_1)
########################################################################################
########################################################################################
# Estimate of the parameters
experiment <- read.table("./wild_pp1_mrna.txt", header=TRUE)
rw <- dim(experiment)[1]
names <- array("", rw)
for (i in 1:rw)
{
names[i] <- "pp1_mrna"
}
names
observed_data_part_1 <- data.frame(name = names,
time = experiment[,1], val = experiment[,2])
observed_data_part_1
ode_solutions_part_1
Cost_function <- function (pars) {
out <- ode_solutions_part_1
cost <- modCost(model = out, obs = observed_data_part_1, y =
"val")
cost
}
Cost_function(parms)
# Fit the model to the observed data
Fit <- modFit(f = Cost_function, p = parms_part_1)
Fit
# Summary of the fit
summary(Fit)
# Model coefficients
coef(Fit)
# Deviance of the fit
deviance(Fit)
