Updated copyright to 2013 for R/ demo/ models/passing and src/ folders, and also...
[openmx:openmx.git] / R / MxExpectationNormal.R
1 #
2 #   Copyright 2007-2013 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 = "MxExpectationNormal",
17         representation = representation(
18                 covariance = "MxCharOrNumber",
19                 means = "MxCharOrNumber",
20                 definitionVars = "list",
21                 thresholds = "MxCharOrNumber",
22                 dims = "character",
23                 dataColumns = "numeric",
24                 thresholdColumns = "numeric",
25                 thresholdLevels = "numeric",
26                 threshnames = "character"),
27         contains = "MxBaseExpectation")
28
29 setMethod("initialize", "MxExpectationNormal",
30         function(.Object, covariance, means, dims, thresholds, threshnames,
31                 data = as.integer(NA), definitionVars = list(), name = 'expectation') {
32                 .Object@name <- name
33                 .Object@covariance <- covariance
34                 .Object@means <- means
35                 .Object@data <- data
36                 .Object@definitionVars <- definitionVars
37                 .Object@thresholds <- thresholds
38                 .Object@dims <- dims
39                 .Object@threshnames <- threshnames
40                 return(.Object)
41         }
42 )
43
44 setMethod("genericExpFunNamespace", signature("MxExpectationNormal"), 
45         function(.Object, modelname, namespace) {
46                 .Object@name <- imxIdentifier(modelname, .Object@name)
47                 .Object@covariance <- imxConvertIdentifier(.Object@covariance, 
48                         modelname, namespace)
49                 .Object@means <- imxConvertIdentifier(.Object@means, 
50                         modelname, namespace)
51                 .Object@data <- imxConvertIdentifier(.Object@data, 
52                         modelname, namespace)
53                 .Object@thresholds <- sapply(.Object@thresholds,
54                         imxConvertIdentifier, modelname, namespace)
55                 return(.Object)
56 })
57
58 setMethod("genericExpDependencies", signature("MxExpectationNormal"),
59         function(.Object, dependencies) {
60         sources <- c(.Object@covariance, .Object@means, .Object@thresholds)
61         sources <- sources[!is.na(sources)]
62         dependencies <- imxAddDependency(sources, .Object@name, dependencies)
63         return(dependencies)
64 })
65
66 setMethod("genericExpAddEntities", "MxExpectationNormal",
67         function(.Object, job, flatJob, labelsData) {
68                 precision <- "Function precision"
69                 if(!single.na(.Object@thresholds)) {
70                         if (is.null(job@options[[precision]])) {
71                                 job <- mxOption(job, precision, 1e-9)
72                         }
73                 }
74                 return(job)
75         }
76 )
77
78 setMethod("genericExpConvertEntities", "MxExpectationNormal",
79         function(.Object, flatModel, namespace, labelsData) {
80                 flatModel <- updateExpectationDimnames(.Object, flatModel, labelsData)
81                 flatModel <- updateThresholdDimnames(.Object, flatModel, labelsData)
82                 return(flatModel)
83         }
84 )
85
86 setMethod("genericExpRename", signature("MxExpectationNormal"),
87         function(.Object, oldname, newname) {
88                 .Object@means <- renameReference(.Object@means, oldname, newname)
89                 .Object@covariance <- renameReference(.Object@covariance, oldname, newname)
90                 .Object@data <- renameReference(.Object@data, oldname, newname)
91                 .Object@thresholds <- sapply(.Object@thresholds, renameReference, oldname, newname)             
92                 return(.Object)
93 })
94
95 verifyExpectedObservedNames <- function(data, covName, flatModel, modelname, objectiveName) {
96         covariance <- flatModel[[covName]]
97         if (is(covariance, "MxMatrix") && !identical(dim(covariance), dim(data))) {
98                 msg <- paste("The dimensions for the expected covariance matrix",
99                         "and the observed covariance matrix",
100                         "in the", objectiveName, "expectation function in model",
101                         omxQuotes(modelname), "are not identical.")
102                 stop(msg, call. = FALSE)
103         }
104         if (!identical(dimnames(covariance), dimnames(data))) {
105                 msg <- paste("The dimnames for the expected covariance matrix",
106                         "and the observed covariance matrix",
107                         "in the", objectiveName, "expectation function in model",
108                         omxQuotes(modelname), "are not identical.")
109                 stop(msg, call. = FALSE)                
110         }
111 }
112
113 verifyMeans <- function(meansName, mxDataObject, flatModel, modelname) {
114         means <- flatModel[[meansName]]
115         if (!is.null(means)) {
116                 if(any(is.na(mxDataObject@means))) {
117                         msg <- paste("In model", omxQuotes(modelname),
118                                 "the Normal expectation function contains an expected means",
119                                 "vector but the model is missing some data",
120                                 "for the observed means.")
121                         stop(msg, call. = FALSE)
122                 }
123                 meanDimnames <- dimnames(means)         
124         }
125 }
126
127 setMethod("genericExpFunConvert", "MxExpectationNormal", 
128         function(.Object, flatModel, model, labelsData, defVars, dependencies) {
129                 modelname <- imxReverseIdentifier(model, .Object@name)[[1]]
130                 name <- .Object@name
131                 if(is.na(.Object@data)) {
132                         msg <- paste("The normal expectation function",
133                                 "does not have a dataset associated with it in model",
134                                 omxQuotes(modelname))
135                         stop(msg, call.=FALSE)
136                 }
137                 mxDataObject <- flatModel@datasets[[.Object@data]]
138                 if (mxDataObject@type != "raw") {
139                         verifyExpectedObservedNames(mxDataObject@observed, .Object@covariance, flatModel, modelname, "Normal")
140                         verifyMeans(.Object@means, mxDataObject, flatModel, modelname)
141                 }
142                 verifyObservedNames(mxDataObject@observed, mxDataObject@means, mxDataObject@type, flatModel, modelname, "Normal")
143                 checkNumericData(mxDataObject)
144                 checkNumberOrdinalColumns(mxDataObject)
145                 meansName <- .Object@means
146                 covName <- .Object@covariance
147                 dataName <- .Object@data
148                 threshName <- .Object@thresholds
149                 covariance <- flatModel[[covName]]
150                 covNames <- colnames(covariance)
151                 .Object@definitionVars <- imxFilterDefinitionVariables(defVars, .Object@data)
152                 .Object@means <- imxLocateIndex(flatModel, .Object@means, name)
153                 .Object@covariance <- imxLocateIndex(flatModel, .Object@covariance, name)
154                 .Object@data <- imxLocateIndex(flatModel, .Object@data, name)
155                 verifyExpectedNames(covName, meansName, flatModel, modelname, "normal")
156                 .Object@dataColumns <- generateDataColumns(flatModel, covNames, dataName)
157                 verifyThresholds(flatModel, model, labelsData, dataName, covNames, threshName)
158                 .Object@thresholds <- imxLocateIndex(flatModel, threshName, name)
159                 retval <- generateThresholdColumns(flatModel, model, labelsData, covNames, dataName, threshName)
160                 .Object@thresholdColumns <- retval[[1]] 
161                 .Object@thresholdLevels <- retval[[2]]
162                 if (length(mxDataObject@observed) == 0) {
163                         .Object@data <- as.integer(NA)
164                 }
165                 if (single.na(.Object@dims)) {
166                         .Object@dims <- covNames
167                 }
168                 return(.Object)
169 })
170
171 verifyExpectedNames <- function(covName, meansName, flatModel, modelname, expectationName) {
172         if (is.na(meansName)) {
173                 means <- NA
174         } else {
175                 means <- flatModel[[meansName]]
176         }
177         covariance <- flatModel[[covName]]
178         covariance <- dimnames(covariance)
179         if (is.null(covariance)) {
180                         msg <- paste("The expected covariance matrix associated",
181                                 "with the", expectationName, "expectation in model", 
182                                 omxQuotes(modelname), "does not contain dimnames.")
183                         stop(msg, call. = FALSE)        
184         }
185         covRows <- covariance[[1]]
186         covCols <- covariance[[2]]      
187         if (is.null(covRows) || is.null(covCols) ||
188                 (length(covRows) != length(covCols)) || !all(covRows == covCols)) {
189                         msg <- paste("The expected covariance matrix associated",
190                                 "with the", expectationName, "expectation in model", 
191                                 omxQuotes(modelname), "does not contain identical",
192                                 "row and column dimnames.")
193                         stop(msg, call.=FALSE)
194         }
195         if (!isS4(means) && is.na(means)) return()
196         meanDimnames <- dimnames(means)
197         if (is.null(meanDimnames)) {
198                         msg <- paste("The expected means matrix associated",
199                                 "with the", expectationName, "expectation function in model", 
200                                 omxQuotes(modelname), "does not contain dimnames.")
201                         stop(msg, call.=FALSE)  
202         }
203         meanRows <- meanDimnames[[1]]
204         meanCols <- meanDimnames[[2]]
205         if (!is.null(meanRows) && length(meanRows) > 1) {
206                 msg <- paste("The expected means matrix associated",
207                         "with the", expectationName, "expectation in model", 
208                         omxQuotes(modelname), "is not a 1 x N matrix.")
209                         stop(msg, call.=FALSE)
210         }
211         if ((length(covCols) != length(meanCols)) || !all(covCols == meanCols)) {
212                         msg <- paste("The expected covariance and expected",
213                                 "means matrices associated",
214                                 "with the", expectationName, "expectation function in model", 
215                                 omxQuotes(modelname), "do not contain identical",
216                                 "dimnames.")
217                         stop(msg, call.=FALSE)
218         }
219 }
220
221 verifyObservedNames <- function(data, means, type, flatModel, modelname, expectationName) {
222         dataNames <- dimnames(data)
223         if(is.null(dataNames)) {
224                 msg <- paste("The observed data associated with the",
225                         expectationName, "expectation function in model",
226                         omxQuotes(modelname), "does not contain dimnames.")
227                 stop(msg, call. = FALSE)
228         }
229         if (type == "cov" || type == "cor") {
230                 if (length(dataNames) < 2 ||
231                         is.null(dataNames[[1]]) || is.null(dataNames[[2]]) || 
232                         !identical(dataNames[[1]], dataNames[[2]])) {
233                                 msg <- paste("The dataset associated with the", expectationName,
234                                         "expectation function in model", omxQuotes(modelname),
235                     "does not contain identical row and column non-NULL dimnames.")
236                         stop(msg, call. = FALSE)
237                 }
238                 if (!single.na(means) && is.null(dimnames(means))) {
239                         msg <- paste("In model", omxQuotes(modelname), 
240                                 ", the observed means vector does not contain column names.",
241                                 "Use the name() function to assign names to the means vector.")
242                         stop(msg, call. = FALSE)
243                 }
244                 if (!single.na(means) && !identical(dataNames[[1]], dimnames(means)[[2]])) {
245                         msg <- paste("The observed covariance or correlation matrix associated with the", expectationName,
246                                 "expectation function in model", omxQuotes(modelname),
247                                 "does not contain identical dimnames to the observed means vector.")
248                         stop(msg, call. = FALSE)
249                 }
250         } else if ((type == "raw") && (length(dataNames) < 2 || is.null(dataNames[[2]]))) {
251                 msg <- paste("The dataset associated with the", expectationName,
252                                 "expectation function in model", omxQuotes(modelname),
253                                 "does not contain column names (use dimnames).")
254                 stop(msg, call. = FALSE)
255         }
256 }
257
258 generateDataColumns <- function(flatModel, covNames, dataName) {
259         retval <- c()
260         dataColumnNames <- dimnames(flatModel@datasets[[dataName]]@observed)[[2]]
261         for(i in 1:length(covNames)) {
262                 targetName <- covNames[[i]]
263                 index <- match(targetName, dataColumnNames)
264                 if(is.na(index)) {
265                         msg <- paste("The column name", omxQuotes(targetName),
266                                 "in the expected covariance matrix",
267                                 "of the expectation function in model",
268                                 omxQuotes(flatModel@name),
269                                 "cannot be found in the column names of the data.")
270                         stop(msg, call. = FALSE)
271                 }
272                 retval[[i]] <- index - 1
273         }
274         return(retval)
275 }
276
277
278 updateThresholdDimnames <- function(flatExpectation, flatModel, labelsData) {
279         threshName <- flatExpectation@thresholds
280         if (is.na(threshName)) {
281                 return(flatModel)
282         }
283         thresholds <- flatModel[[threshName]]
284         if (is.null(thresholds)) {
285                 modelname <- getModelName(flatExpectation)
286                 stop(paste("Unknown thresholds name", 
287                         omxQuotes(simplifyName(threshName, modelname)),
288                         "detected in the expectation function",
289                         "of model", omxQuotes(modelname)), call. = FALSE)
290         }
291         dims <- flatExpectation@threshnames
292         if (!is.null(colnames(thresholds)) && !single.na(dims) && 
293                 !identical(colnames(thresholds), dims)) {
294                 modelname <- getModelName(flatExpectation)
295                 msg <- paste("The thresholds matrix associated",
296                 "with the expectation function in model", 
297                 omxQuotes(modelname), "contains column names and",
298                 "the expectation function has specified non-identical threshnames.")
299                 stop(msg, call.=FALSE)      
300         }
301         if (is.null(colnames(thresholds)) && !single.na(dims)) {
302                 tuple <- evaluateMxObject(threshName, flatModel, labelsData, new.env(parent = emptyenv()))
303                 threshMatrix <- tuple[[1]]
304                 if (ncol(threshMatrix) != length(dims)) {
305                         modelname <- getModelName(flatExpectation)
306                         msg <- paste("The thresholds matrix associated",
307                         "with the expectation function in model", 
308                         omxQuotes(modelname), "is not of the same length as the 'threshnames'",
309                         "argument provided by the expectation function. The 'threshnames' argument is",
310                         "of length", length(dims), "and the expected covariance matrix",
311                         "has", ncol(threshMatrix), "columns.")
312                         stop(msg, call.=FALSE)      
313                 }
314                 dimnames(flatModel[[threshName]]) <- list(NULL, dims)
315         }
316         return(flatModel)
317 }
318
319 updateExpectationDimnames <- function(flatExpectation, flatModel,
320                 labelsData, unsafe = FALSE) {
321         covName <- flatExpectation@covariance
322         meansName <- flatExpectation@means
323         if (is.na(meansName)) {
324                 means <- NA
325         } else {
326                 means <- flatModel[[meansName]]
327         }
328         covariance <- flatModel[[covName]]
329         if (is.null(covariance)) {
330                 modelname <- getModelName(flatExpectation)
331                 stop(paste("Unknown expected covariance name", 
332                         omxQuotes(simplifyName(covName, modelname)),
333                         "detected in the expectation function",
334                         "of model", omxQuotes(modelname)), call. = FALSE)
335         }
336         if (is.null(means)) {
337                 modelname <- getModelName(flatExpectation)
338                 stop(paste("Unknown expected means name", 
339                         omxQuotes(simplifyName(meansName, modelname)),
340                         "detected in the expectation function",
341                         "of model", omxQuotes(modelname)), call. = FALSE)
342         }
343         dims <- flatExpectation@dims
344         if (!is.null(dimnames(covariance)) && !single.na(dims) && 
345                 !identical(dimnames(covariance), list(dims, dims))) {
346                 modelname <- getModelName(flatExpectation)
347                 msg <- paste("The expected covariance matrix associated",
348                         "with the expectation function in model", 
349                         omxQuotes(modelname), "contains dimnames: ", 
350             paste(toString(dimnames(covariance)), ".", sep = ""),
351                         "The expectation function has specified dimnames:", 
352                         paste(toString(dims), ".", sep =""))
353                 stop(msg, call.=FALSE)          
354         }
355         if (is.null(dimnames(covariance)) && !single.na(dims)) {
356                 if (!unsafe) {
357                         tuple <- evaluateMxObject(covName, flatModel, labelsData, new.env(parent = emptyenv()))
358                         covMatrix <- tuple[[1]]
359                         if (nrow(covMatrix) != ncol(covMatrix)) {
360                                 modelname <- getModelName(flatExpectation)
361                                 msg <- paste("The expected covariance matrix associated",
362                                         "with the expectation function in model", 
363                                         omxQuotes(modelname), "is not a square matrix.")
364                                 stop(msg, call.=FALSE)          
365                         }
366                         if (nrow(covMatrix) != length(dims)) {
367                                 modelname <- getModelName(flatExpectation)
368                                 msg <- paste("The expected covariance matrix associated",
369                                         "with the expectation function in model", 
370                                         omxQuotes(modelname), "is not of the same length as the 'dimnames'",
371                                         "argument provided by the expectation function. The 'dimnames' argument is",
372                                         "of length", length(dims), "and the expected covariance matrix",
373                                         "has", nrow(covMatrix), "rows and columns.")
374                                 stop(msg, call.=FALSE)          
375                         }
376                 }
377                 dimnames(flatModel[[covName]]) <- list(dims, dims)
378         }
379
380         if (!isS4(means) && is.na(means)) {
381                 return(flatModel)
382         }
383
384         if (!is.null(dimnames(means)) && !single.na(dims) &&
385                 !identical(dimnames(means), list(NULL, dims))) {
386                 modelname <- getModelName(flatExpectation)
387                 msg <- paste("The expected means matrix associated",
388                         "with the expectation function in model", 
389                         omxQuotes(modelname), "contains dimnames: ", 
390             paste(toString(dimnames(means)), ".", sep = ""),
391                         "The expectation function has specified dimnames:", 
392                         paste(toString(dims), ".", sep =""))
393                 stop(msg, call.=FALSE)  
394         }
395         if (is.null(dimnames(means)) && !single.na(dims)) {
396                 if (!unsafe) {
397                         tuple <- evaluateMxObject(meansName, flatModel, labelsData, new.env(parent = emptyenv()))
398                         meansMatrix <- tuple[[1]]
399                         if (nrow(meansMatrix) != 1) {
400                                 modelname <- getModelName(flatExpectation)
401                                 msg <- paste("The expected means vector associated",
402                                         "with the expectation function in model", 
403                                         omxQuotes(modelname), "is not a 1 x n matrix.",
404                                         "It has dimensions", nrow(meansMatrix), "x", 
405                                         paste(ncol(meansMatrix), '.', sep=''))
406                                 stop(msg, call.=FALSE)          
407                         }
408                         if (ncol(meansMatrix) != length(dims)) {
409                                 modelname <- getModelName(flatExpectation)
410                                 msg <- paste("The expected means vector associated",
411                                         "with the expectation function in model", 
412                                         omxQuotes(modelname), "is not of the same length as the 'dimnames'",
413                                         "argument provided by the expectation function. The 'dimnames' argument is",
414                                         "of length", length(dims), "and the expected means vector",
415                                         "has", ncol(meansMatrix), "columns.")
416                                 stop(msg, call.=FALSE)
417                         }
418                 }
419                 dimnames(flatModel[[meansName]]) <- list(NULL, dims)
420         }
421         return(flatModel)
422 }
423
424
425
426 mxExpectationNormal <- function(covariance, means = NA, 
427         dimnames = NA, thresholds = NA, threshnames = dimnames) {
428         if (missing(covariance) || typeof(covariance) != "character") {
429                 stop("'covariance' argument is not a string (the name of the expected covariance matrix)")
430         }
431         if (!(single.na(means) || typeof(means) == "character")) {
432                 stop("Means argument is not a string (the name of the expected means matrix)")
433         }
434         if (is.na(means)) means <- as.integer(NA)
435         if (single.na(thresholds)) thresholds <- as.character(NA)
436         if (single.na(dimnames)) dimnames <- as.character(NA)
437         if (single.na(threshnames)) threshnames <- as.character(NA)
438         if (!is.vector(dimnames) || typeof(dimnames) != 'character') {
439                 stop("'dimnames' argument is not a character vector")
440         }
441         if (!is.vector(threshnames) || typeof(threshnames) != 'character') {
442                 stop("'threshnames' argument is not a character vector")
443         }
444         if (length(thresholds) != 1) {
445                 stop("'thresholds' argument must be a single matrix or algebra name")
446         }
447         if (length(dimnames) == 0) {
448                 stop("'dimnames' argument cannot be an empty vector")
449         }
450         if (length(threshnames) == 0) {
451                 stop("'threshnames' argument cannot be an empty vector")
452         }
453         if (length(dimnames) > 1 && any(is.na(dimnames))) {
454                 stop("NA values are not allowed for 'dimnames' vector")
455         }
456         if (length(threshnames) > 1 && any(is.na(threshnames))) {
457                 stop("NA values are not allowed for 'threshnames' vector")
458         }
459         return(new("MxExpectationNormal", covariance, means, dimnames, thresholds, threshnames))
460 }
461
462 displayMxExpectationNormal <- function(expectation) {
463         cat("MxExpectationNormal", omxQuotes(expectation@name), '\n')
464         cat("@covariance :", omxQuotes(expectation@covariance), '\n')
465         cat("@means :", omxQuotes(expectation@means), '\n')
466         if (single.na(expectation@dims)) {
467                 cat("@dims : NA \n")
468         } else {
469                 cat("@dims :", omxQuotes(expectation@dims), '\n')
470         }
471         if (single.na(expectation@thresholds)) {
472                 cat("@thresholds : NA \n")
473         } else {
474                 cat("@thresholds :", omxQuotes(expectation@thresholds), '\n')
475         }
476         if (single.na(expectation@threshnames)) {
477                 cat("@threshnames : NA \n")
478         } else {
479                 cat("@threshnames :", omxQuotes(expectation@threshnames), '\n')
480         }
481         cat("@info$likelihoods: ", length(expectation@info$likelihoods) > 0, '\n')
482         invisible(expectation)
483 }
484
485
486 setMethod("print", "MxExpectationNormal", function(x,...) { 
487         displayMxExpectationNormal(x) 
488 })
489
490 setMethod("show", "MxExpectationNormal", function(object) { 
491         displayMxExpectationNormal(object) 
492 })