Implemented global variable substitution in mxAlgebra statements.
[openmx:openmx.git] / R / MxModelFunctions.R
1 #
2 #   Copyright 2007-2009 The OpenMx Project
3 #
4 #   Licensed under the Apache License, Version 2.0 (the "License");
5 #   you may not use this file except in compliance with the License.
6 #   You may obtain a copy of the License at
7
8 #        http://www.apache.org/licenses/LICENSE-2.0
9
10 #   Unless required by applicable law or agreed to in writing, software
11 #   distributed under the License is distributed on an "AS IS" BASIS,
12 #   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 #   See the License for the specific language governing permissions and
14 #   limitations under the License.
15
16
17 generateMatrixList <- function(mxModel) {
18         matvalues <- lapply(mxModel@matrices, generateMatrixValuesHelper)
19         matnames  <- names(mxModel@matrices)
20         names(matvalues) <- matnames
21         references <- generateMatrixReferences(mxModel)
22         retval <- mapply(function(x,y) { c(list(x), y) }, 
23                         matvalues, references, SIMPLIFY = FALSE)
24         return(retval)
25 }
26
27 generateSimpleMatrixList <- function(mxModel) {
28         retval <- generateMatrixList(mxModel)
29         retval <- lapply(retval, function(x) { 
30                 c(list(as.matrix(x[[1]])), x[-1]) }) 
31         return(retval)
32 }
33
34 generateAlgebraList <- function(mxModel) {
35         mNames <- names(mxModel@matrices)
36         aNames <- names(mxModel@algebras)
37         oNames <- names(mxModel@objectives)
38         retval <- lapply(mxModel@algebras, generateAlgebraHelper, 
39                 mNames, append(aNames, oNames))
40     return(retval)
41 }
42
43 generateParameterList <- function(flatModel) {
44         result <- list()
45         if (length(flatModel@matrices) == 0) {
46                 return(result)
47         }
48         for(i in 1:length(flatModel@matrices)) {
49                 matrix <- flatModel@matrices[[i]]
50                 result <- generateParameterListHelper(
51                         matrix, result, i - 1L)
52         }
53         return(result)
54 }
55
56 generateDefinitionList <- function(flatModel) {
57         result <- list()
58         defLocations <- generateDefinitionLocations(flatModel@datasets)
59         if (length(flatModel@matrices) == 0) {
60                 return(result)
61         }
62         for(i in 1:length(flatModel@matrices)) {
63                 result <- generateDefinitionListHelper(
64                         flatModel@matrices[[i]], 
65                         result, defLocations, i - 1L)
66         }       
67         return(result)
68 }
69
70 generateValueList <- function(mxModel, mList, pList) {
71         mList <- lapply(mList, function(x) { x[[1]] })
72         retval <- vector()
73         if (length(pList) == 0) {
74                 return(retval)
75         }
76         for(i in 1:length(pList)) {
77                 parameter <- pList[[i]]
78                 parameter <- parameter[3:length(parameter)] # Remove (min, max) bounds
79                 if (length(parameter) > 1) {
80                         values <- sapply(parameter, generateValueHelper, mList)
81                         if (!all(values == values[[1]])) {
82                                 warning(paste('Parameter',i,'has multiple start values.',
83                                         'Selecting', values[[1]]))
84                         }
85                         retval[i] <- values[[1]]
86                 } else {
87                         retval[i] <- generateValueHelper(parameter[[1]], mList)
88                 }
89     }
90         return(retval)  
91 }
92
93 generateValueHelper <- function(triple, mList) {
94         mat <- triple[1] + 1
95         row <- triple[2] + 1
96         col <- triple[3] + 1
97         return(mList[[mat]][row,col])
98 }
99
100 getObjectiveIndex <- function(flatModel) {
101         objective <- flatModel@objective
102         if(is.null(objective)) {
103                 return(NULL)
104         } else {
105                 return(omxLocateIndex(flatModel, objective@name, flatModel@name))
106         }
107 }
108
109 updateModelValues <- function(model, flatModel, pList, values) {
110         if(length(pList) != length(values)) {
111                 stop(paste("This model has", length(pList), 
112                         "parameters, but you have given me", length(values),
113                         "values"))
114         }
115         if (length(pList) == 0) {
116                 return(model)
117         }
118         for(i in 1:length(pList)) {
119                 parameters <- pList[[i]]
120                 parameters <- parameters[3:length(parameters)] # Remove (min, max) bounds
121                 model <- updateModelValuesHelper(parameters, values[[i]], flatModel, model)
122     }
123         return(model)
124 }
125
126 updateModelValuesHelper <- function(triples, value, flatModel, model) {
127         for(i in 1:length(triples)) {
128                 triple <- triples[[i]]
129                 mat <- triple[1] + 1
130                 row <- triple[2] + 1
131                 col <- triple[3] + 1
132                 name <- flatModel@matrices[[mat]]@name
133                 model[[name]]@values[row,col] <- value
134         }
135         return(model)
136 }
137
138 removeTail <- function(lst, tailSize) {
139     newEnd <- length(lst) - tailSize
140     if (newEnd == 0) {
141         return(list())
142     } else {
143         return(lst[1 : newEnd])
144     }
145 }
146
147 updateModelMatrices <- function(model, flatModel, values) {
148     flatModel@matrices <- removeTail(flatModel@matrices, length(flatModel@constMatrices))
149     flatModel@matrices <- removeTail(flatModel@matrices, length(flatModel@freeMatrices))
150     flatModel@matrices <- removeTail(flatModel@matrices, length(flatModel@outsideMatrices))
151     values <- removeTail(values, length(flatModel@constMatrices))    
152     values <- removeTail(values, length(flatModel@freeMatrices))
153     values <- removeTail(values, length(flatModel@outsideMatrices))    
154         mList <- names(flatModel@matrices)
155         if (length(mList) != length(values)) {
156                 stop(paste("This model has", length(mList), 
157                         "matrices, but you have given me", length(values),
158                         "values"))
159         }
160         if (length(mList) == 0) {
161                 return(model)
162         }       
163         model <- updateModelMatricesHelper(mList, values, model)
164         return(model)
165 }
166
167
168 updateModelAlgebras <- function(model, flatModel, values) {
169         aNames <- names(flatModel@algebras)
170         oNames <- names(flatModel@objectives)
171         aList <- append(aNames, oNames)
172         if(length(aList) != length(values)) {
173                 stop(paste("This model has", length(aList), 
174                         "algebras, but you have given me", length(values),
175                         "values"))
176         }
177         if (length(aList) == 0) {
178                 return(model)
179         }       
180         model <- updateModelAlgebrasHelper(aList, values, model)
181         return(model)
182 }
183
184 updateModelMatricesHelper <- function(mList, values, model) {
185         for(i in 1:length(mList)) {
186                 name <- mList[[i]]
187                 model[[name]]@values <- values[[i]]
188         }
189         return(model)
190 }
191
192 updateModelAlgebrasHelper <- function(aList, values, model) {
193         for(i in 1:length(aList)) {
194                 name <- aList[[i]]
195                 candidate <- model[[name]]
196                 if (!is.null(candidate) && (length(values[[i]]) > 0) 
197                         && !is.nan(values[[i]]) && 
198                         (is(candidate,"MxAlgebra") || (is(candidate,"MxObjective")))) {
199                         model[[name]]@result <- as.matrix(values[[i]])
200                         if (is(candidate, "MxAlgebra")) {
201                                 dimnames(model[[name]]@result) <- dimnames(model[[name]])
202                         }
203                 }
204         }
205         return(model)
206 }
207
208 omxLocateIndex <- function(model, name, referant) {
209         if (is.na(name)) { return(as.integer(name)) }
210         mNames <- names(model@matrices)
211         aNames <- names(model@algebras)
212         oNames <- names(model@objectives)
213         dNames <- names(model@datasets)         
214         matrixNumber <- match(name, mNames)
215         algebraNumber <- match(name, append(aNames, oNames))
216         dataNumber <- match(name, dNames)
217         if (is.na(matrixNumber) && is.na(algebraNumber) && is.na(dataNumber)) {
218                 msg <- paste("The reference", omxQuotes(name),
219                         "does not exist.  It is used by the named entity",
220                         omxQuotes(referant),".")
221                 stop(msg, call.=FALSE)
222         } else if (!is.na(matrixNumber)) {
223                 return(- matrixNumber)
224         } else if (!is.na(dataNumber)) {
225                 return(dataNumber - 1L)
226         } else {
227                 return(algebraNumber - 1L)
228         }
229 }
230
231 omxCheckMatrices <- function(model) {
232         matrices <- model@matrices
233         lapply(matrices, omxVerifyMatrix)
234         if (length(model@submodels) > 0) {
235                 for(i in 1:length(model@submodels)) {
236                         submodel <- model@submodels[[i]]
237                         if(submodel@independent == FALSE) {
238                                 omxCheckMatrices(submodel)
239                         }
240                 }
241         }
242 }