Rework specification of free param groups per Tim's suggestion
[openmx:openmx.git] / R / MxCompute.R
1 #
2 #   Copyright 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 = "MxBaseCompute", 
17          representation = representation(
18            id = "integer",
19            "VIRTUAL"),
20          contains = "MxBaseNamed")
21
22 setClassUnion("MxCompute", c("NULL", "MxBaseCompute"))
23
24 setGeneric("convertForBackend",
25         function(.Object, flatModel, model) {
26                 return(standardGeneric("convertForBackend"))
27         })
28
29 setGeneric("assignId",
30         function(.Object, id) {
31                 return(standardGeneric("assignId"))
32         })
33
34 setMethod("assignId", signature("MxBaseCompute"),
35         function(.Object, id) {
36                 .Object@id <- id
37                 .Object
38         })
39
40 setGeneric("getFreeVarGroup",
41         function(.Object) {
42                 return(standardGeneric("getFreeVarGroup"))
43         })
44
45 setMethod("getFreeVarGroup", signature("MxBaseCompute"),
46         function(.Object) {
47                 list()
48         })
49
50 #----------------------------------------------------
51
52 setClass(Class = "MxComputeOperation",
53          contains = "MxBaseCompute",
54          representation = representation(
55            free.set = "MxOptionalChar"))
56
57 setMethod("qualifyNames", signature("MxComputeOperation"),
58         function(.Object, modelname, namespace) {
59                 .Object@name <- imxIdentifier(modelname, .Object@name)
60                 .Object@free.set <- imxConvertIdentifier(.Object@free.set, modelname, namespace)
61                 .Object
62         })
63
64 setMethod("getFreeVarGroup", signature("MxComputeOperation"),
65         function(.Object) {
66                 if (length(.Object@free.set)) {
67                         list(.Object@id, .Object@free.set)
68                 } else {
69                         list()
70                 }
71         })
72
73 setMethod("convertForBackend", signature("MxComputeOperation"),
74         function(.Object, flatModel, model) {
75                 name <- .Object@name
76                 .Object
77         })
78
79 #----------------------------------------------------
80
81 setClass(Class = "MxComputeOnce",
82          contains = "MxComputeOperation",
83          representation = representation(
84            what = "MxCharOrNumber",
85            context = "character",
86            gradient = "logical",
87            hessian = "logical"))
88
89 setMethod("qualifyNames", signature("MxComputeOnce"),
90         function(.Object, modelname, namespace) {
91                 .Object <- callNextMethod();
92                 .Object@what <- imxConvertIdentifier(.Object@what, modelname, namespace)
93                 .Object
94         })
95
96 setMethod("convertForBackend", signature("MxComputeOnce"),
97         function(.Object, flatModel, model) {
98                 .Object <- callNextMethod();
99                 name <- .Object@name
100                 if (any(!is.integer(.Object@what))) {
101                         expNum <- match(.Object@what, names(flatModel@expectations))
102                         algNum <- match(.Object@what, append(names(flatModel@algebras),
103                                                              names(flatModel@fitfunctions)))
104                         if (any(is.na(expNum)) && any(is.na(algNum))) {
105                                 stop("Can only apply MxComputeOnce to MxAlgebra or MxExpectation")
106                         }
107                         if (!any(is.na(expNum))) {
108                                         # Usually negative numbers indicate matrices; not here
109                                 .Object@what <- - expNum
110                         } else {
111                                 if (any(algNum > length(flatModel@algebras)) && length(algNum) > 1) {
112                                         stop("MxComputeOnce cannot evaluate more than 1 fit function")
113                                 }
114                                 .Object@what <- algNum - 1L
115                         }
116                 }
117                 if (length(.Object@what) == 0) warning("MxComputeOnce with nothing will have no effect")
118                 .Object
119         })
120
121 setMethod("initialize", "MxComputeOnce",
122           function(.Object, what, free.set, context, gradient, hessian) {
123                   .Object@name <- 'compute'
124                   .Object@what <- what
125                   .Object@free.set <- free.set
126                   .Object@context <- context
127                   .Object@gradient <- gradient
128                   .Object@hessian <- hessian
129                   .Object
130           })
131
132 mxComputeOnce <- function(what, free.set=NULL, context=character(0), gradient=FALSE, hessian=FALSE) {
133         new("MxComputeOnce", what, free.set, context, gradient, hessian)
134 }
135
136 #----------------------------------------------------
137
138 setClass(Class = "MxComputeGradientDescent",
139          contains = "MxComputeOperation",
140          representation = representation(
141            fitfunction = "MxCharOrNumber",
142            engine = "character"))
143
144 setMethod("qualifyNames", signature("MxComputeGradientDescent"),
145         function(.Object, modelname, namespace) {
146                 .Object <- callNextMethod();
147                 .Object@fitfunction <- imxConvertIdentifier(.Object@fitfunction, modelname, namespace)
148                 .Object
149         })
150
151 setMethod("convertForBackend", signature("MxComputeGradientDescent"),
152         function(.Object, flatModel, model) {
153                 .Object <- callNextMethod();
154                 name <- .Object@name
155                 if (is.character(.Object@fitfunction)) {
156                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
157                 }
158                 .Object
159         })
160
161 setMethod("initialize", "MxComputeGradientDescent",
162           function(.Object, free.set, engine, fit) {
163                   .Object@name <- 'compute'
164                   .Object@free.set <- free.set
165                   .Object@fitfunction <- fit
166                   .Object@engine <- engine
167                   .Object
168           })
169
170 mxComputeGradientDescent <- function(type, free.set=NULL,
171                                      engine=NULL, fitfunction='fitfunction') {
172 # What to do with 'type'?
173 #       if (length(type) != 1) stop("Specific 1 compute type")
174
175         if (is.null(engine)) engine <- as.character(NA)
176
177         new("MxComputeGradientDescent", free.set, engine, fitfunction)
178 }
179
180 #----------------------------------------------------
181
182 setClass(Class = "MxComputeNewtonRaphson",
183          contains = "MxComputeOperation",
184          representation = representation(
185            fitfunction = "MxCharOrNumber",
186            maxIter = "integer",
187            tolerance = "numeric"))
188
189 setMethod("qualifyNames", signature("MxComputeNewtonRaphson"),
190         function(.Object, modelname, namespace) {
191                 .Object <- callNextMethod();
192                 .Object@fitfunction <- imxConvertIdentifier(.Object@fitfunction, modelname, namespace)
193                 .Object
194         })
195
196 setMethod("convertForBackend", signature("MxComputeNewtonRaphson"),
197         function(.Object, flatModel, model) {
198                 .Object <- callNextMethod();
199                 name <- .Object@name
200                 if (is.character(.Object@fitfunction)) {
201                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
202                 }
203                 .Object
204         })
205
206 setMethod("initialize", "MxComputeNewtonRaphson",
207           function(.Object, free.set, fit, maxIter, tolerance) {
208                   .Object@name <- 'compute'
209                   .Object@free.set <- free.set
210                   .Object@fitfunction <- fit
211                   .Object@maxIter <- maxIter
212                   .Object@tolerance <- tolerance
213                   .Object
214           })
215
216 mxComputeNewtonRaphson <- function(type, free.set=NULL,
217                                    fitfunction='fitfunction', maxIter = 500L, tolerance=1e-7) {
218
219         new("MxComputeNewtonRaphson", free.set, fitfunction, maxIter, tolerance)
220 }
221
222 #----------------------------------------------------
223
224 setClass(Class = "MxComputeSteps",
225          contains = "MxBaseCompute",
226          representation = representation(
227            steps = "list"))
228
229 setMethod("getFreeVarGroup", signature("MxComputeSteps"),
230         function(.Object) {
231                 result <- list()
232                 for (step in .Object@steps) {
233                         got <- getFreeVarGroup(step)
234                         if (length(got)) result <- append(result, got)
235                 }
236                 result
237         })
238
239 setMethod("assignId", signature("MxComputeSteps"),
240         function(.Object, id) {
241                 steps <- .Object@steps
242                 for (sx in 1:length(steps)) {
243                         steps[[sx]] <- assignId(steps[[sx]], id)
244                         id <- steps[[sx]]@id + 1L
245                 }
246                 .Object@steps <- steps
247                 .Object@id <- id
248                 .Object
249         })
250
251 setMethod("qualifyNames", signature("MxComputeSteps"),
252         function(.Object, modelname, namespace) {
253                 .Object@name <- imxIdentifier(modelname, .Object@name)
254                 .Object@steps <- lapply(.Object@steps, function (c) qualifyNames(c, modelname, namespace))
255                 .Object
256         })
257
258 setMethod("convertForBackend", signature("MxComputeSteps"),
259         function(.Object, flatModel, model) {
260                 .Object@steps <- lapply(.Object@steps, function (c) convertForBackend(c, flatModel, model))
261                 .Object
262         })
263
264 #----------------------------------------------------
265
266 setClass(Class = "MxComputeIterate",
267          contains = "MxComputeSteps",
268          representation = representation(
269            maxIter = "integer",
270            tolerance = "numeric",
271            verbose = "logical"))
272
273 setMethod("initialize", "MxComputeIterate",
274           function(.Object, steps, maxIter, tolerance, verbose) {
275                   .Object@name <- 'compute'
276                   .Object@steps <- steps
277                   .Object@maxIter <- maxIter
278                   .Object@tolerance <- tolerance
279                   .Object@verbose <- verbose
280                   .Object
281           })
282
283 mxComputeIterate <- function(steps, maxIter=500L, tolerance=1e-4, verbose=FALSE) {
284         new("MxComputeIterate", steps=steps, maxIter=maxIter, tolerance=tolerance, verbose)
285 }
286
287 displayMxComputeIterate <- function(opt) {
288         cat(class(opt), omxQuotes(opt@name), '\n')
289         cat("@tolerance :", omxQuotes(opt@tolerance), '\n')
290         cat("@maxIter :", omxQuotes(opt@maxIter), '\n')
291         for (step in 1:length(opt@steps)) {
292                 cat("[[", step, "]] :", class(opt@steps[[step]]), '\n')
293         }
294         invisible(opt)
295 }
296
297 setMethod("print", "MxComputeIterate", function(x, ...) displayMxComputeIterate(x))
298 setMethod("show",  "MxComputeIterate", function(object) displayMxComputeIterate(object))
299
300 #----------------------------------------------------
301
302 setClass(Class = "MxComputeEstimatedHessian",
303          contains = "MxComputeOperation",
304          representation = representation(
305            fitfunction = "MxCharOrNumber",
306            se = "logical"))
307
308 setMethod("qualifyNames", signature("MxComputeEstimatedHessian"),
309         function(.Object, modelname, namespace) {
310                 .Object <- callNextMethod();
311                 .Object@fitfunction <- imxConvertIdentifier(.Object@fitfunction, modelname, namespace)
312                 .Object
313         })
314
315 setMethod("convertForBackend", signature("MxComputeEstimatedHessian"),
316         function(.Object, flatModel, model) {
317                 .Object <- callNextMethod();
318                 name <- .Object@name
319                 if (is.character(.Object@fitfunction)) {
320                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
321                 }
322                 .Object
323         })
324
325 setMethod("initialize", "MxComputeEstimatedHessian",
326           function(.Object, free.set, fit, want.se) {
327                   .Object@name <- 'compute'
328                   .Object@free.set <- free.set
329                   .Object@fitfunction <- fit
330                   .Object@se <- want.se
331                   .Object
332           })
333
334 mxComputeEstimatedHessian <- function(free.set=NULL, fitfunction='fitfunction', want.se=TRUE) {
335         new("MxComputeEstimatedHessian", free.set, fitfunction, want.se)
336 }
337
338 #----------------------------------------------------
339
340 setClass(Class = "MxComputeSequence",
341          contains = "MxComputeSteps")
342
343 setMethod("initialize", "MxComputeSequence",
344           function(.Object, steps) {
345                   .Object@name <- 'compute'
346                   .Object@steps <- steps
347                   .Object
348           })
349
350 mxComputeSequence <- function(steps) {
351         new("MxComputeSequence", steps=steps)
352 }
353
354 displayMxComputeSequence <- function(opt) {
355         cat(class(opt), omxQuotes(opt@name), '\n')
356         for (step in 1:length(opt@steps)) {
357                 cat("[[", step, "]] :", class(opt@steps[[step]]), '\n')
358         }
359         invisible(opt)
360 }
361
362 setMethod("print", "MxComputeSequence", function(x, ...) displayMxComputeSequence(x))
363 setMethod("show",  "MxComputeSequence", function(object) displayMxComputeSequence(object))
364
365 #----------------------------------------------------
366
367 displayMxComputeOperation <- function(opt) {
368         cat(class(opt), omxQuotes(opt@name), '\n')
369         cat("@id :", opt@id, '\n')
370         cat("@free.set :", omxQuotes(opt@free.set), '\n')
371         invisible(opt)
372 }
373
374 setMethod("print", "MxComputeOperation", function(x, ...) displayMxComputeOperation(x))
375 setMethod("show",  "MxComputeOperation", function(object) displayMxComputeOperation(object))
376
377 displayMxComputeGradientDescent <- function(opt) {
378         cat("@type :", omxQuotes(opt@type), '\n')
379         cat("@engine :", omxQuotes(opt@engine), '\n')
380         cat("@fitfunction :", omxQuotes(opt@fitfunction), '\n')
381         invisible(opt)
382 }
383
384 setMethod("print", "MxComputeGradientDescent",
385           function(x, ...) { callNextMethod(); displayMxComputeGradientDescent(x) })
386 setMethod("show",  "MxComputeGradientDescent",
387           function(object) { callNextMethod(); displayMxComputeGradientDescent(object) })
388
389 convertComputes <- function(flatModel, model) {
390         retval <- lapply(flatModel@computes, function(opt) {
391                 convertForBackend(opt, flatModel, model)
392         })
393         retval
394 }