Added the data field to an MxModel.
[openmx:openmx.git] / R / MxModel.R
1 setClass(Class = "MxModel",
2         representation = representation(
3                 matrices = "list",
4                 algebras = "list",
5                 data = "data.frame"))
6                 
7 setMethod("initialize", "MxModel",
8         function(.Object, matrices = list(), algebras = list(), data = data.frame()) {
9                 .Object@matrices = matrices
10                 .Object@algebras = algebras
11                 .Object@data = data
12                 return(.Object)
13         }
14 )
15
16 setMethod("[[", "MxModel",
17         function(x, i, j, ..., drop = FALSE) {
18                 first <- x@matrices[[i]]
19                 second <- x@algebras[[i]]
20                 if (is.null(first)) {
21                         return(second)
22                 } else {
23                         return(first)
24                 }       
25         }
26 )
27
28 setReplaceMethod("[[", "MxModel",
29         function(x, i, j, value) {
30                 if (is(value,"MxMatrix")) {
31                         if (!is.null(x@algebras[[i]])) {
32                                 stop(paste(i, "is already an MxAlgebra object"))
33                         }
34                         x@matrices[[i]] <- value
35                 } else if (is(value,"MxAlgebra")) {
36                         if (!is.null(x@matrices[[i]])) {
37                                 stop(paste(i, "is already an MxMatrix object"))
38                         }
39                         x@algebras[[i]] <- value                
40                 } else {
41                         stop(paste("Unknown type of value", value))
42                 }
43                 return(x)
44         }
45 )
46
47 omxGenerateMatrixList <- function(mxModel) {
48         return(lapply(mxModel@matrices, generateMatrixListHelper))
49 }
50
51 omxGenerateSimpleMatrixList <- function(mxModel) {
52         retval <- lapply(mxModel@matrices, generateMatrixListHelper)
53         return(lapply(retval, as.matrix))
54 }
55
56 omxGenerateAlgebraList <- function(mxModel) {
57     mList <- omxGenerateMatrixList(mxModel)
58     retval <- lapply(mxModel@algebras, generateAlgebraHelper, names(mList))
59     return(retval)
60 }
61
62 omxGenerateParameterList <- function(mxModel) {
63         result <- list()
64         for(i in 1:length(mxModel@matrices)) {
65                 result <- generaterParameterListHelper(mxModel@matrices[[i]], result, i - 1)
66         }       
67         return(result)
68 }
69
70 omxGenerateValueList <- function(mxModel) {
71         mList <- omxGenerateMatrixList(mxModel)
72         pList <- omxGenerateParameterList(mxModel)
73         retval <- vector()
74         for(i in 1:length(pList)) {
75                 parameter <- pList[[i]]
76                 if (length(parameter) > 1) {
77                         values <- sapply(parameter, generateValueHelper, mList)
78                         if (!all(values == values[[1]])) {
79                                 warning(paste('Parameter',i,'has multiple start values.',
80                                         'Selecting', values[[1]]))
81                         }
82                         retval[i] <- values[[1]]
83                 } else {
84                         retval[i] <- generateValueHelper(parameter[[1]], mList)
85                 }
86     }
87         return(retval)  
88 }
89
90 options(stringsAsFactors = FALSE)
91 data(omxSymbolTable)
92 options(stringsAsFactors = TRUE)
93
94 generateAlgebraHelper <- function(algebra, lnames) {
95         retval <- algebra@formula
96         numbers <- as.list(as.double(-1 : (-length(lnames))))
97         names(numbers) <- lnames
98         retval <- eval(substitute(substitute(e, numbers), list(e = retval)))
99         retval <- substituteOperators(retval)
100         return(retval)
101 }
102
103 substituteOperators <- function(algebra) {
104         if ((length(algebra) > 1) && (!is.numeric(algebra[[1]]))) {
105                 names <- omxSymbolTable["R.name"] == as.character(algebra[[1]])
106         variableSymbols <- omxSymbolTable["Number.of.arguments"] == -1
107                 result <- omxSymbolTable[names & variableSymbols, "Num"]
108                 if (length(result) > 1) {
109                                 stop(paste("Ambiguous function with name", algebra[[1]],
110                                         "and", (length(algebra) - 1), "arguments"))
111                 } else if(length(result) == 1) {
112                         head <- as.double(result[[1]])
113                         tail <- lapply(algebra[-1], substituteOperators)
114                         return(list(head, tail))
115                 } else {
116                         length <- omxSymbolTable["Number.of.arguments"] == (length(algebra) - 1)
117                         result <- omxSymbolTable[names & length, "Num"]
118                         if (length(result) == 0) {
119                                 stop(paste("Could not find function with name", algebra[[1]],
120                                         "and", (length(algebra) - 1), "arguments"))
121                         } else if (length(result) > 1) {
122                                 stop(paste("Ambiguous function with name", algebra[[1]],
123                                         "and", (length(algebra) - 1), "arguments"))
124                         } else {
125                                 head <- as.double(result[[1]])
126                             tail <- lapply(algebra[-1], substituteOperators)
127                                 result <- append(tail, head, after=0)
128                                 return(result)
129                         }
130         }
131         }
132         return(algebra)
133 }
134
135 generateValueHelper <- function(triple, mList) {
136         mat <- triple[1] + 1
137         row <- triple[2]
138         col <- triple[3]
139         return(mList[[mat]][row,col])
140 }
141
142 omxUpdateModelValues <- function(mxModel, values) {
143         pList <- omxGenerateParameterList(mxModel)
144         if(length(pList) != length(values)) {
145                 stop(paste("This model has", length(pList), 
146                         "parameters, but you have given me", length(values),
147                         "values"))
148         }
149         for(i in 1:length(pList)) {
150                 mxModel <- updateModelValueHelper(pList[[i]], values[[i]], mxModel)
151     }
152         return(mxModel)
153 }
154
155 updateModelValueHelper <- function(triples, value, mxModel) {
156         for(i in 1:length(triples)) {
157                 triple <- triples[[i]]
158                 mat <- triple[1] + 1
159                 row <- triple[2]
160                 col <- triple[3]
161                 mxModel@matrices[[mat]]@values[row,col] <- value                        
162         }
163         return(mxModel)
164 }