Fixed a bug where could not be accessed from a model object, and adjusted one test...
[openmx:openmx.git] / R / MxModel.R
1 #
2 #   Copyright 2007-2015 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 setClass(Class = "MxModel",
17         representation = representation(
18                 name = "character",
19                 matrices = "list",
20                 algebras = "list",
21                 constraints = "list",
22                 intervals = "list",
23                 latentVars = "MxCharOrList",
24                 manifestVars = "MxCharOrList",
25                 data = "MxData",
26                 submodels = "list",
27                 expectation = "MxExpectation",
28                 fitfunction = "MxFitFunction",
29                 compute = "MxCompute",
30                 independent = "logical",
31                 options = "list",
32                 output = "list",
33                 runstate = "list",
34                 .forcesequential = "logical",
35                 .newobjects = "logical",
36                 .resetdata = "logical",
37                 .wasRun = "logical",
38                 .modifiedSinceRun = "logical"
39 ))
40
41 imxModelTypes[['default']] <- "MxModel"
42
43 setMethod("initialize", "MxModel",
44         function(.Object, name = character()) {
45                 .Object@name <- name
46                 .Object@latentVars <- character()
47                 .Object@manifestVars <- character()
48                 .Object@matrices <- list()
49                 .Object@algebras <- list()
50                 .Object@constraints <- list()
51                 .Object@data <- NULL
52                 .Object@submodels <- list()
53                 .Object@expectation <- NULL
54                 .Object@fitfunction <- NULL
55                 .Object@compute <- NULL
56                 .Object@independent <- FALSE
57                 .Object@options <- list()
58                 .Object@output <- list()
59                 .Object@runstate <- list()
60                 .Object@.newobjects <- FALSE
61                 .Object@.resetdata <- FALSE
62                 .Object@.wasRun <- FALSE
63                 .Object@.modifiedSinceRun <- FALSE
64                 .Object <- imxInitModel(.Object)
65                 return(.Object)
66         }
67 )
68
69 ##' imxInitModel
70 ##'
71 ##' This is an internal function exported for those people who know
72 ##' what they are doing.
73 ##'
74 ##' @param model model
75 ##' @aliases
76 ##' imxInitModel,MxModel-method
77 ##' imxInitModel,MxRAMModel-method
78 ##' imxInitModel,MxLISRELModel-method
79 setGeneric("imxInitModel", function(model) {
80         return(standardGeneric("imxInitModel")) } )
81
82 ##' imxModelBuilder
83 ##'
84 ##' This is an internal function exported for those people who know
85 ##' what they are doing.
86 ##' 
87 ##' TODO: It probably makes sense to split this into separate
88 ##' methods. For example, modelAddVariables and modelRemoveVariables
89 ##' could be their own methods. This would reduce some cut&paste
90 ##' duplication.
91 ##'
92 ##' @param model model
93 ##' @param lst lst
94 ##' @param name name
95 ##' @param manifestVars manifestVars
96 ##' @param latentVars latentVars
97 ##' @param submodels submodels
98 ##' @param remove remove
99 ##' @param independent independent
100 ##' @aliases
101 ##' imxModelBuilder,MxLISRELModel-method
102 ##' imxModelBuilder,MxModel-method
103 ##' imxModelBuilder,MxRAMModel-method
104 setGeneric("imxModelBuilder", function(model, lst, name, 
105         manifestVars, latentVars, submodels, remove, independent) {
106         return(standardGeneric("imxModelBuilder")) } )
107
108 ##' imxTypeName
109 ##'
110 ##' This is an internal function exported for those people who know
111 ##' what they are doing.
112 ##'
113 ##' @param model model
114 ##' @aliases
115 ##' imxTypeName,MxLISRELModel-method
116 ##' imxTypeName,MxModel-method
117 ##' imxTypeName,MxRAMModel-method
118 setGeneric("imxTypeName", function(model) { 
119         return(standardGeneric("imxTypeName")) 
120 })
121
122 ##' imxVerifyModel
123 ##' 
124 ##' This is an internal function exported for those people who know
125 ##' what they are doing.
126 ##'
127 ##' @param model model
128 ##' @aliases
129 ##' imxVerifyModel,MxLISRELModel-method
130 ##' imxVerifyModel,MxModel-method
131 ##' imxVerifyModel,MxRAMModel-method
132 setGeneric("imxVerifyModel", function(model) {
133     return(standardGeneric("imxVerifyModel"))
134 })
135
136 # End declaration of generics
137
138 generateParentNames <- function(model) {
139         retval <- generateLocalNames(model)
140         if (length(model@submodels) > 0) {
141                 retval <- union(retval, names(model@submodels))
142                 childNames <- unlist(lapply(model@submodels, generateChildNames))
143                 retval <- union(retval, childNames)
144         }
145         return(retval)
146 }
147
148 generateChildNames <- function(model) {
149         retval <- generateLocalNames(model)     
150         if (!is.null(retval)) {
151                 retval <- paste(model@name, retval, sep = ".")
152         }
153         if (length(model@submodels) > 0) {
154                 retval <- union(retval, names(model@submodels))
155                 childNames <- unlist(lapply(model@submodels, generateChildNames))
156                 retval <- union(retval, childNames)
157         }
158         return(retval)
159 }
160
161 generateLocalNames <- function(model) {
162         matrices <- names(model@matrices)
163         algebras <- names(model@algebras)
164         constraints <- names(model@constraints)
165         retval <- union(matrices, algebras)
166         retval <- union(retval, constraints)
167         if (!is.null(model@fitfunction)) {
168                 retval <- union(retval, model@fitfunction@name)
169         }
170         if (!is.null(model@expectation)) {
171                 retval <- union(retval, model@expectation@name)
172         }
173         if (!is.null(model@data)) {
174                 retval <- union(retval, model@data@name)
175         }
176         return(retval)
177 }
178
179 setMethod("[[", "MxModel",
180         function(x, i, j, ..., drop = FALSE) {
181                 return(imxExtractMethod(x, i))
182         }
183 )
184
185 setReplaceMethod("[[", "MxModel",
186         function(x, i, j, value) {
187                 return(imxReplaceMethod(x, i, value))
188         }
189 )
190
191 # These are slots that are intended to be directly viewable by the user.
192 # Included separately so that they are the same between the $ and names() operators.
193 publicMxModelSlots <- c("name", "matrices", "algebras", "constraints", "data", "submodels", "output", "compute", "options", "intervals")
194 visibleMxModelSlots <- c("name", "options", "compute", "output", "intervals")
195
196 setMethod("$", "MxModel",
197         function(x, name) {
198         result <- imxExtractMethod(x, name)
199         if(name %in% publicMxModelSlots) {
200             result <- imxExtractSlot(x, name)
201         }
202                 return(result)
203         }
204 )
205
206 setMethod("names", "MxModel",
207         function(x) {
208                 submodels <- names(x@submodels)
209                 locals <- generateLocalNames(x)
210                 slots <- publicMxModelSlots
211                 output <- c(submodels, locals, slots)
212                 output <- gsub("(\\w+\\W+.*)", "'\\1'", output)
213                 return(output)
214         }
215 )
216
217 .DollarNames.MxModel <- function(x, pattern) {   # .DollarNames is an S3 Generic
218                 submodels <- names(x@submodels)
219                 locals <- generateLocalNames(x)
220                 slots <- imxGetSlotDisplayNames(x, slotList=visibleMxModelSlots)
221                 output <- c(submodels, locals, slots)
222                 output <- gsub("(\\w+\\W+.*)", "'\\1'", output)
223                 return(grep(pattern, output, value=TRUE))
224         }
225
226 setReplaceMethod("$", "MxModel",
227         function(x, name, value) {
228                 if(name == "output") {
229                         stop("You cannot directly set the output of a model.  Use mxRun() if you want output.")
230                 } 
231                 if(name == "name") {
232                         stop("You cannot directly set the name of a model.  To rename the model, use model<-mxModel(model, name=\"NewName\").")
233                 }
234                 if(name %in% c("matrices", "algebras", "submodels")) {
235                         stop(paste("You cannot directly set the", name,
236                                    "of a model.  To set objects in the model, use the mxModel() function."))
237                 }
238                 return(imxReplaceMethod(x, name, value))
239         }
240 )
241
242
243 ##' imxSameType
244 ##' 
245 ##' This is an internal function exported for those people who know
246 ##' what they are doing.
247 ##'
248 ##' @param a a
249 ##' @param b b
250 imxSameType <- function(a, b) {
251         return( (is(a, "MxModel") && is(b, "MxModel")) ||
252                         (is(a, "MxMatrix") && is(b, "MxMatrix")) ||
253                         (is(a, "MxAlgebra") && is(b, "MxAlgebra")) ||
254                         (is(a, "MxExpectation") && is(b, "MxExpectation")) ||
255                         (is(a, "MxFitFunction") && is(b, "MxFitFunction")) ||
256                         (is(a, "MxCompute") && is(b, "MxCompute")) ||
257                         (is(a, "MxConstraint") && is(b, "MxConstraint")) ||
258                         (is(a, "MxData") && is(b, "MxData")))
259 }
260
261 mxModel <- function(model = NA, ..., manifestVars = NA, latentVars = NA,
262                     remove = FALSE, independent = NA, type = NA, name = NA) {
263         retval <- firstArgument(model, name)
264         first <- retval[[1]]
265         model <- retval[[2]]
266         name  <- retval[[3]]
267         model <- typeArgument(model, type)
268         lst <- c(first, list(...))
269         lst <- unlist(lst)
270         filter <- sapply(lst, is, "MxModel")
271         submodels <- lst[filter]
272         lst <- lst[!filter]
273         model <- imxModelBuilder(model, lst, name, manifestVars,
274                 latentVars, submodels, remove, independent)
275         return(model)
276 }
277
278 firstArgument <- function(model, name) {
279         first <- NULL
280         defaultType <- imxModelTypes[[getOption("mxDefaultType")]]
281         if (is(model, "MxModel")) {
282         } else {
283                 if (single.na(model)) {
284                 } else if (typeof(model) == "character") {
285                         name <- model
286                 } else if (isS4(model)) {
287                         first <- model
288                 } else {
289                         first <- list(model)
290                 }
291                 if (length(name) > 0 && is.na(name)) {
292                         name <- imxUntitledName()
293                 }
294                 imxVerifyName(name, -1)
295                 model <- new(defaultType, name)
296         }
297         return(list(first, model, name))
298 }
299
300 typeArgument <- function(model, type) {
301         if (!is.na(type)) {
302                 if (is.null(imxModelTypes[[type]])) {
303                         stop(paste("The model type", omxQuotes(type), 
304                                 "is not in the the list of acceptable types:",
305                                 omxQuotes(names(imxModelTypes))), call. = FALSE)
306                 }
307                 typename <- imxModelTypes[[type]]
308                 class(model) <- typename
309                 model <- imxInitModel(model)
310         }
311         return(model)
312 }
313
314 ##' imxGenericModelBuilder
315 ##'
316 ##' This is an internal function exported for those people who know
317 ##' what they are doing.
318 ##'
319 ##' @param model model
320 ##' @param lst lst
321 ##' @param name name
322 ##' @param manifestVars manifestVars
323 ##' @param latentVars latentVars
324 ##' @param submodels submodels
325 ##' @param remove remove
326 ##' @param independent independent
327 imxGenericModelBuilder <- function(model, lst, name, 
328         manifestVars, latentVars, submodels, remove, independent) {
329         model <- nameArgument(model, name)
330         model <- variablesArgument(model, manifestVars, latentVars, submodels, remove)
331         model <- listArgument(model, lst, remove)
332         model <- independentArgument(model, independent)
333         return(model)
334 }
335
336 varsToCharacter <- function(vars, vartype) {
337         if (is.list(vars)) {
338                 varnames <- names(vars)
339                 if (length(varnames) == 0) {
340                         return(as.character(unlist(vars)))      
341                 } else {
342                         result <- pmatch(varnames, imxVariableTypes)
343                         illegal <- which(is.na(result))
344                         if (length(illegal) > 0) {
345                                 if (length(illegal) == 1) {
346                                         ctgMsg <- "category"
347                                 } else {
348                                         ctgMsg <- "categories"
349                                 }
350                                 msg <- paste("In the", vartype, "variables",
351                                         "the", ctgMsg,
352                                         omxQuotes(varnames[illegal]), "did not match",
353                                         "to a valid category or two categories matched",
354                                         "to the same string (see 'imxVariableTypes'",
355                                         "for the list of legal categories)")
356                                 stop(msg, call. = FALSE)
357                         }
358                         varnames <- imxVariableTypes[result]
359                         vars <- lapply(vars, as.character)
360                         names(vars) <- varnames
361                         return(vars)
362                 }
363         } else {
364                 return(as.character(vars))
365         }
366 }
367
368 variablesArgument <- function(model, manifestVars, latentVars, submodels, remove) {
369         if (single.na(manifestVars)) {
370                 manifestVars <- character()
371         }
372         if (single.na(latentVars)) {
373                 latentVars <- character()
374         }
375         if (remove == TRUE) {
376                 model <- modelRemoveVariables(model, latentVars, manifestVars)
377                 if (length(submodels)) for(i in 1:length(submodels)) {
378                         model <- removeSingleNamedEntity(model, submodels[[i]])
379                 }
380         } else {
381                 if (length(manifestVars) + length(latentVars) > 0) {
382                         latentVars <- varsToCharacter(latentVars, "latent")
383                         manifestVars <- varsToCharacter(manifestVars, "manifest")
384                         checkVariables(model, latentVars, manifestVars)
385                         model <- modelAddVariables(model, latentVars, manifestVars)
386                 }
387                 if (length(submodels)) for(i in 1:length(submodels)) {
388                         model <- addSingleNamedEntity(model, submodels[[i]])
389                 }
390         }
391         return(model)
392 }
393
394 listArgument <- function(model, lst, remove) {
395         if(remove == TRUE) {
396                 model <- modelRemoveEntries(model, lst)
397         } else {
398                 model <- modelAddEntries(model, lst)
399         }
400         return(model)
401 }
402
403 independentArgument <- function(model, independent) {
404         if(!is.na(independent)) {
405                 model@independent <- independent
406         }
407         return(model)
408 }
409
410 nameArgument <- function(model, name) {
411         if(!is.na(name)) {
412                 model@name <- name
413         }
414         return(model)
415 }
416
417 checkVariables <- function(model, latentVars, manifestVars, submodels) {
418         latentVars   <- unlist(latentVars, use.names = FALSE)
419         manifestVars <- unlist(manifestVars, use.names = FALSE)
420         modelLatent <- unlist(model@latentVars, use.names = FALSE)
421         modelManifest <- unlist(model@manifestVars, use.names = FALSE)
422         common <- intersect(latentVars, manifestVars)    
423         if (length(common) > 0) {
424                 stop(paste("The following variables cannot",
425                         "be both latent and manifest:",
426                         omxQuotes(common)), call. = FALSE)
427         }
428         common <- intersect(modelLatent, manifestVars)
429         if (length(common) > 0) {
430                 stop(paste("The following variables cannot",
431                         "be both latent and manifest:",
432                         omxQuotes(common)), call. = FALSE)
433         }
434         common <- intersect(modelManifest, latentVars)
435         if (length(common) > 0) {
436                 stop(paste("The following variables cannot",
437                         "be both latent and manifest",
438                         omxQuotes(common)), call. = FALSE)
439         }
440         common <- intersect(modelManifest, manifestVars)
441         if (length(common) > 0) {
442                 stop(paste("The following manifest variables",
443                         "have already been declared",
444                         omxQuotes(common)), call. = FALSE)
445         }
446         common <- intersect(modelLatent, latentVars)
447         if (length(common) > 0) {
448                 stop(paste("The following latent variables",
449                         "have already been declared",
450                         omxQuotes(common)), call. = FALSE)
451         }
452         if (any(is.na(latentVars))) {
453                 stop("NA is not allowed as a latent variable", call. = FALSE)
454         }
455         if (any(is.na(manifestVars))) {
456                 stop("NA is not allowed as a manifest variable", call. = FALSE)
457         }
458         if (length(unique(latentVars)) != length(latentVars)) {
459                 stop(paste("The following variables in the latentVars list are duplicated:", 
460                 omxQuotes(latentVars[duplicated(latentVars)])), call. = FALSE)
461         }
462         if (length(unique(manifestVars)) != length(manifestVars)) {
463                 stop(paste("The following variables in the manifestVars list are duplicated:", 
464                 omxQuotes(manifestVars[duplicated(manifestVars)])), call. = FALSE)
465         }
466 }
467
468 # Begin implementation of generics
469
470 setMethod("imxModelBuilder", "MxModel", imxGenericModelBuilder)
471
472 setMethod("imxInitModel", "MxModel", function(model) { 
473         return(model)
474 })
475
476 setMethod("imxTypeName", "MxModel", function(model) { 
477         return("default")
478 })
479
480 setMethod("imxVerifyModel", "MxModel", function(model) {
481     return(TRUE)
482 })
483
484 # End implementation of generics
485
486 addVariablesHelper <- function(model, vartype, vars) {
487         modelvars <- slot(model, vartype)
488
489         if (length(vars) == 0) {
490                 return(model)
491         } else if (length(modelvars) == 0) {
492                 slot(model, vartype) <- vars
493                 return(model)
494         }
495
496         if (is.list(vars) && !is.list(modelvars)) {
497                 msg <- paste("The", vartype, "variables in",
498                         "the call to mxModel() have been separated",
499                         "into categories, and the existing", vartype,
500                         "variables do not have categories.")
501                 stop(msg, call. = FALSE)
502         } else if (!is.list(vars) && is.list(modelvars)) {
503                 msg <- paste("The", vartype, "variables in",
504                         "the call to mxModel() have not been separated",
505                         "into categories, and the existing", vartype,
506                         "variables do have categories.")
507                 stop(msg, call. = FALSE)
508         }
509
510         if (is.character(vars) && is.character(modelvars)) {
511                 modelvars <- c(modelvars, vars)
512                 slot(model, vartype) <- modelvars
513         } else {
514                 varnames <- names(vars)
515                 offsets <- pmatch(varnames, imxVariableTypes)
516                 for(i in 1:length(vars)) {
517                         currOffset <- offsets[[i]]
518                         currName <- imxVariableTypes[[currOffset]]
519                         modelvars[[currName]] <- c(modelvars[[currName]], vars[[i]])
520                 }
521                 slot(model, vartype) <- modelvars
522         }
523
524         return(model)
525 }
526
527 modelAddVariables <- function(model, latent, manifest) {
528         model <- addVariablesHelper(model, "latentVars", latent)
529         model <- addVariablesHelper(model, "manifestVars", manifest)
530         return(model)
531 }
532
533 modelRemoveVariables <- function(model, latent, manifest) {
534         model@latentVars <- setdiff(model@latentVars, latent)
535         model@manifestVars <- setdiff(model@manifestVars, manifest)
536         return(model)
537 }
538         
539 modelAddEntries <- function(model, entries) {
540         if (length(entries) == 0) {
541                 return(model)
542         }
543         tuple <- modelModifyFilter(model, entries, "add")
544         namedEntities <- tuple[[1]]
545         bounds        <- tuple[[2]]
546         intervals     <- tuple[[3]]
547         intervals     <- expandIntervals(intervals)
548         names(intervals) <- sapply(intervals, slot, "reference")
549         if (length(namedEntities) > 0) for(i in 1:length(namedEntities)) {
550                 model <- addSingleNamedEntity(model, namedEntities[[i]])
551         }
552         model <- modelAddBounds(model, bounds)
553         model <- modelAddIntervals(model, intervals)
554         return(model)
555 }
556
557 modelRemoveEntries <- function(model, entries) {
558         if (length(entries) == 0) {
559                 return(model)
560         }
561         tuple <- modelModifyFilter(model, entries, "remove")
562         namedEntities <- tuple[[1]]
563         bounds        <- tuple[[2]]
564         intervals     <- tuple[[3]]
565         intervals     <- expandIntervals(intervals)     
566         names(intervals) <- sapply(intervals, slot, "reference")
567         if (length(namedEntities) > 0) for(i in 1:length(namedEntities)) {
568                 model <- removeSingleNamedEntity(model, namedEntities[[i]])
569         }
570         model <- modelRemoveBounds(model, bounds)
571         model <- modelRemoveIntervals(model, intervals)
572         return(model)
573 }
574
575 actionCorrespondingPredicate <- c('add' = 'into', 'remove' = 'from')
576
577 modelModifyFilter <- function(model, entries, action) {
578         boundsFilter <- sapply(entries, is, "MxBounds")
579         intervalFilter <- sapply(entries, is, "MxInterval")
580         namedEntityFilter <- sapply(entries, function(x) {"name" %in% slotNames(x)})
581         characterFilter <- sapply(entries, is.character)
582         pathFilter <- sapply(entries, is, "MxPath")
583   thresholdFilter <- sapply(entries, is, "MxThreshold")
584         unknownFilter <- !(boundsFilter | namedEntityFilter | intervalFilter | characterFilter | thresholdFilter)
585         if (any(pathFilter)) {
586                 stop(paste("The model type of model",
587                         omxQuotes(model@name), "does not recognize paths."),
588                         call. = FALSE)
589         }
590         if (any(thresholdFilter)) {
591           stop(paste("The model type of model",
592                      omxQuotes(model@name), "does not recognize thresholds."),
593                call. = FALSE)
594         }
595         if (any(unknownFilter)) {
596                 stop(paste("Cannot", action, "the following item(s)", 
597                         actionCorrespondingPredicate[[action]], "the model:", 
598                         omxQuotes(sapply(entries[unknownFilter], deparse))), call. = FALSE)
599         }
600         if (any(namedEntityFilter) && action == 'remove') {
601                 stop(paste("Cannot use named entities when remove = TRUE.",
602                         "Instead give the name of the entity when removing it.",
603                         "See http://openmx.psyc.virginia.edu/wiki/mxmodel-help#Remove_an_object_from_a_model"))
604         }
605         if (any(characterFilter) && action == 'add') {
606                 stop(paste("I don't know what to do with the following strings",
607                         omxQuotes(entries[characterFilter]),
608                         "that have been passed into the function:",
609                         deparse(width.cutoff = 400L, imxLocateFunction("mxModel"))), call. = FALSE)
610         }
611         if (identical(action, 'add')) {
612                 return(list(entries[namedEntityFilter], entries[boundsFilter], entries[intervalFilter]))
613         } else if (identical(action, 'remove')) {
614                 return(list(entries[characterFilter], entries[boundsFilter], entries[intervalFilter]))
615         } else {
616                 stop(paste("Internal error, unidentified action:", omxQuotes(action)))
617         }
618 }
619
620 addSingleNamedEntity <- function(model, entity) {
621         if (!nzchar(entity@name)) {
622                 stop(paste("Entity",
623                         omxQuotes(class(entity)), "in model",
624                         omxQuotes(model@name), "needs a name"), call. = FALSE)
625         }
626         if (model@name == entity@name) {
627                 stop(paste("You cannot insert an entity named",
628                         omxQuotes(entity@name), "into a model named",
629                         omxQuotes(model@name)), call. = FALSE)
630         }
631         model[[entity@name]] <- entity
632         return(model)
633 }
634
635 removeSingleNamedEntity <- function(model, name) {
636         model[[name]] <- NULL
637         return(model)
638 }
639
640 setMethod("imxVerifyModel", "MxModel",
641     function(model) {
642         if (length(model@submodels) > 0) {
643                 return(all(sapply(model@submodels, imxVerifyModel)))
644         }
645         return(TRUE)
646     }
647 )