changes to doc files
[openmx:openmx.git] / docs / source / NewBeginnersGuide.rst
1 Beginners Guide to OpenMx
2 =========================
3
4 This document is a basic introduction to the OpenMx library.  It assumes that the reader has installed the R statistical programming language [http://www.r-project.org/] and the OpenMx library for R [http://openmx.psyc.virginia.edu].  Detailed introductions to R can be found on the internet.  We recommend [http://faculty.washington.edu/tlumley/Rcourse] for a short course but Google search for 'introduction to R' provides many options.
5
6 The OpenMx scripts for the examples in this guide are available in the following files:
7
8 * http://openmx.psyc.virginia.edu/repoview/1/trunk/demo/OneFactorPathDemo.R
9 * http://openmx.psyc.virginia.edu/repoview/1/trunk/demo/OneFactorMatrixDemo.R
10
11 * http://www.vipbg.vcu.edu/OpenMx/html/NewBeginnersGuide.R
12
13
14 Basic Introduction 
15 ------------------
16
17 The main function of OpenMx is to provide a way to design a statistical model which can be used to test a hypothesis with a set of data.  The way to do this is to use one of the functions written in the R language as part of the OpenMx package.  The function to create an Mx model is (you guessed it) ``mxModel()``.  Note that i) R is case sensitive and ii) the OpenMx package must be loaded before any of the OpenMx functions are used (this only has to be done once in an R session).
18
19 .. cssclass:: input
20 ..
21
22    OpenMx Code
23    
24    .. code-block:: r
25        
26         require(OpenMx)
27
28 We will start by building a model with the ``mxModel()`` function, which takes a number of arguments.  We will explain each of these gradually below.  Usually we save the result of applying the ``mxModel()`` function as an R object, here called *myModel*.  This R object has a special class named MxModel.
29
30 .. cssclass:: input
31 ..
32
33    OpenMx Code
34    
35    .. code-block:: r
36        
37         myModel <- mxModel() 
38
39 We might read the above line 'myModel gets mxModel left paren right paren'. In this case we have provided no arguments in the parentheses, but R has still created an empty MxModel object *myModel*. Obviously this model needs arguments to do anything useful, but just to get the idea of the process of 'Model Building - Model Fitting' here is how we would go about fitting this model.  The *myModel* object becomes the argument of the OpenMx function ``mxRun()`` to fit the model.  The result is stored in a new R object, *myModelRun* of the same class as the previous object but updated with the results from the model fitting.
40
41 .. cssclass:: input
42 ..
43
44    OpenMx Code
45    
46    .. code-block:: r
47        
48         myModelRun <- mxRun(myModel) 
49
50
51 A model typically contains data, free and/or fixed parameters and some function to be applied to the data according to the model.  Models can be build in a variety of ways, due to the flexibility of OpenMx which it inherited from classic **Mx** - for those of you who familiar with that software - and further expanded.
52
53 One possibility for building models is to create all the elements as separate R objects, which are then combined by calling them up in the ``mxModel()`` statement.  We will refer to this approach as the **piecewise style**.
54
55 A second type is to start by creating an ``mxModel()`` with one element, and then taking this model as the basis for the next model where you add another element to it and so on.  This also provides a good way to make sure that each element is syntactically correct before adding anything new.  We will call this approach the **iterative/recursive/stepwise style**.
56
57 Another approach more closely resembles the traditional classic **Mx** approach, where you specify all the elements of the model at once, all as arguments of the ``mxModel()``, separated by comma's.  While this approach is often more compact and works well for scripts that run successfully, it is not the most recommended approach for debugging a new script.  We will refer to this approach as the **classic style**.
58
59
60 A First mxModel
61 ----------------
62
63 To introduce the model fitting process in OpenMx, we will present the basics of several OpenMx functions which can be used to write a simple model and view its contents.
64
65 Matrix Creation
66 ^^^^^^^^^^^^^^^
67
68 Although ``mxModel()`` can have a range of arguments, we will start with the most simple one.  Models are fitted to data which must be in numeric format (for continuous data) or factor format (for ordinal data).  Here we consider continuous data.  Numbers (data/parameter estimates) are typically put into matrices, except for fixed constants.  The function created to put numbers into matrices is (unsurprisingly) ``mxMatrix()``.  Here we start with a basic matrix call and make use of only some of its possible arguments. All arguments are separated by comma's. To make it clear and explicit, we will include the names of the arguments, although that is optional if the.arguments are included in the default order.
69
70 .. cssclass:: input
71 ..
72
73    OpenMx Code
74    
75    .. code-block:: r
76        
77         myAmatrix <- mxMatrix(type="Full", nrow=1, ncol=1, values=4myModel, name="Amatrix")
78     
79     
80 The above call to the ``mxMatrix()`` function has five arguments.  The ``type`` and ``name`` arguments are alphanumeric and therefore their values are in quotes.  The ``nrows``, ``ncols`` and ``values`` arguments are numeric, and refer respectively to the number of rows, the number of columns of the matrix and the value for the (in this case only one) element of the matrix.
81
82 Matrix Contents
83 ^^^^^^^^^^^^^^^
84
85 Once you have run/executed this statement in R, a new R object has been created, namely *myAmatrix*.  When you view its contents, you'll notice it has a special class of object, made by OpenMx, called an MxMatrix object.  This object has a number of attributes, all of which are listed when you call up the object.  
86
87     ..  code-block:: r
88
89         > myAmatrix
90         FullMatrix 'Amatrix' 
91         
92         $labels: No labels assigned.
93         
94         $values
95           [,1]
96         [1,]    4
97         
98         $free: No free parameters.
99         
100         $lbound: No lower bounds assigned.
101         
102         $ubound: No upper bounds assigned.
103
104
105 Most of these attributes start with the ``$`` symbol.  The contents of a particular attribute can be displayed by typing the name of the R object followed by the ``$`` symbol and the name of the attribute, for example here we're displaying the values of the matrix *myAmatrix*
106    
107    .. code-block:: r
108    
109         > myAmatrix$values
110                [,1]
111           [1,]    4
112
113
114 Note that the attribute ``name`` is part of the header of the output but is not displayed as an ``$`` attribute.  However, it does exist as one and can be seen by typing
115    
116    .. code-block:: r
117    
118         > myAmatrix$name
119         [1] "Amatrix"
120
121
122 Wait a minute, this is confusing.  The matrix has a name, here "Amatrix", and the R object to represent the matrix has a name, here "myAmatrix".  Remember that when you call up *myAmatrix* you get the contents of the entire MxMatrix R object.  When you call up "Amatrix", you get 
123
124     .. code-block:: r
125
126         Error: object 'Amatrix' not found   
127     
128 unless you had previously created another R object with that same name.  Why do we need two names?  The matrix name (here, "Amatrix") is used within OpenMx when performing an operation on this matrix using algebra (see below) or manipulating/using the matrix in any way within a model.  When you want to manipulate/use/view the matrix outside of OpenMx, or build a model by building each of the elements as R objects in the 'piecewise' approach, you use the R object name (here, *myAmatrix*).  Let's clarify this with an example.  
129
130 Model Creation
131 ^^^^^^^^^^^^^^
132
133 First, we will build a model *myModel1* with just one matrix.  Obviously that is not very useful but it does serve to introduce the sequence of creating a model and running it.  
134
135 .. cssclass:: input
136 ..
137
138    OpenMx Code
139    
140    .. code-block:: r
141
142         myModel1 <- mxModel( 
143                         mxMatrix(type="Full", nrow=1, ncol=1, values=4, name="Amatrix") 
144                     )
145                     
146 Model Execution
147 ^^^^^^^^^^^^^^^^
148
149 The ``mxRun()`` function will run a model through the optimizer.  The return value of this function is an identical MxModel object, with all the free parameters - in case there are any - in the elements of the matrices of the model assigned to their final values.                    
150                     
151 .. cssclass:: input
152 ..
153
154    OpenMx Code
155    
156    .. code-block:: r
157    
158         myModel1Run <- mxRun(myModel1)
159
160 Model Contents
161 ^^^^^^^^^^^^^^
162
163 Note that we have saved the result of applying ``mxRun()`` to *myModel1* into a new R object, called *myModel1Run* which is of the same class as *myModel1* but with values updated after fitting the model.  Note that the MxModel is automatically given a name 'untitled2' as we did not specify a ``name`` argument for the ``mxModel()`` function.
164    
165    .. code-block:: r
166     
167         >     myModel1Run
168         MxModel 'untitled2' 
169         type : default 
170         $matrices : 'Amatrix' 
171         $algebras :  
172         $constraints :  
173         $intervals :  
174         $latentVars : none
175         $manifestVars : none
176         $data : NULL
177         $submodels :  
178         $expectation : NULL 
179         $fitfunction : NULL 
180         $compute : NULL 
181         $independent : FALSE 
182         $options :  
183         $output : TRUE
184
185 As you can see from viewing the contents of the new object, the current model only uses two of the arguments, namely ``$matrices`` and ``$output``.  Given the matrix was specified within the mxModel, we can explore its arguments by extending the level of detail as follows.
186
187    .. code-block:: r
188
189         > myModel1Run$matrices
190         $Amatrix
191         FullMatrix 'Amatrix' 
192     
193         $labels: No labels assigned.
194     
195         $values
196              [,1]
197         [1,]    4
198     
199         $free: No free parameters.
200     
201         $lbound: No lower bounds assigned.
202     
203         $ubound: No upper bounds assigned.
204     
205 This lists all the matrices within the MxModel *myModel1Run*.  In the current case there is only one.  If we want to display just a specific argument of that matrix, we first add a dollar sign, followed by the name of the matrix, and an ``$`` symbol prior to the required argument.  Thus arguments within an object are preceded by the ``$`` symbol; specific elements of the same argument type are preceded by the ``$`` sign.
206
207     .. code-block:: r
208
209         > myModel1run$matrices$Amatrix$values
210               [,1]
211          [1,]    4
212
213 It is also possible to omit the ``$matrices`` part and use the more succinct ``myModel1Run$Amatrix$values``.
214
215 Similarly, we can inspect the output which also includes the matrices in ``$matrices``, but only displays the values.  Furthermore, the output will list algebras (``$algebras``), model expectations (``$expectations``), status of optimization (``$status``), number of evaluations (``$evaluations``), openmx version (``$mxVersion``), and a series of time measures of which the CPU time might be most useful (``$cpuTime``).
216
217     .. code-block:: r
218
219         > myModel1Run$output
220         $matrices
221         $matrices$untitled2.Amatrix
222              [,1]
223         [1,]    4
224
225         $algebras
226         list()
227
228         $expectations
229         list()
230
231         $status
232         $status$code
233         [1] -1
234
235         $status$status
236         [1] 0
237
238
239         $iterations
240         [1] 0
241
242         $evaluations
243         [1] 1 1
244
245         $mxVersion
246         [1] "999.0.0-3297"
247
248         $frontendTime
249         Time difference of 0.05656791 secs
250
251         $backendTime
252         Time difference of 0.003615141 secs
253
254         $independentTime
255         Time difference of 3.385544e-05 secs
256
257         $wallTime
258         Time difference of 0.0602169 secs
259
260         $timestamp
261         [1] "2014-04-10 09:53:37 EDT"
262
263         $cpuTime
264         Time difference of 0.0602169 secs
265         
266
267 Alternative 
268 ^^^^^^^^^^^
269
270 Now let's go back to the model *myModel1* for a minute.  We specified the matrix "Amatrix" within the model.  Given we had previously saved the "Amatrix" in the *myAmatrix* object, we could have just used the R object as the argument of the model as follows.  Here we're adding one additional element to the ``MxModel()`` object, namely the ``name`` argument
271
272 .. cssclass:: input
273 ..
274
275    OpenMx Code
276    
277    .. code-block:: r
278
279         myModel2 <- mxModel(myAmatrix, name="model2")
280     
281         myModel2Run <- mxRun(myModel2)
282
283 You can verify for yourself that the contents of *myModel2* is identical to that of *myModel1*, and the same applies to *myModel1Run* and *myModel2Run*, and as a result to the matrix contained in the model.  The value of the matrix element is still 4, both in the original model and the fitted model, as we did not manipulate the matrix in any way.  We refer to this alternative style of coding as 'iterative'.
284
285 Algebra Creation
286 ^^^^^^^^^^^^^^^^
287
288 Now, let's take it one step further and use OpenMx to evaluate some matrix algebra.  It will come as a bit of a shock to learn that the OpenMx function to specify an algebra is called ``mxAlgebra()``.  Its main argument is the ``expression``, in other words the matrix algebra formula you want to evaluate.  In this case, we're simply adding 1 to the value of the matrix element, providing a name for the matrix "Bmatrix" and then save the new matrix as *myBmatrix*.  Note that the matrix we are manipulating is the "Amatrix", the name given to the matrix within OpenMx.
289
290 .. cssclass:: input
291 ..
292
293    OpenMx Code
294    
295    .. code-block:: r
296
297         myBmatrix <- mxAlgebra(expression=Amatrix+1, name="Bmatrix")
298     
299 Algebra Contents
300 ^^^^^^^^^^^^^^^^
301
302 We can view the contents of this new matrix. Notice that the result has not yet computed, as we have not run the model yet.
303     
304    .. code-block:: r
305     
306         > myBmatrix
307         mxAlgebra 'Bmatrix' 
308         $formula:  Amatrix + 1 
309         $result: (not yet computed) <0 x 0 matrix>
310         dimnames: NULL
311
312
313 Built Model
314 ^^^^^^^^^^^
315
316 Now we can combine the two statements - one defining the matrix, and the other defining the algebra - in one model, simply by separating them by a comma, and run it to see the result of the operation.
317
318 .. cssclass:: input
319 ..
320
321    OpenMx Code
322    
323    .. code-block:: r
324
325         myModel3 <- mxModel(myAmatrix,myBmatrix, name="model3")
326     
327         myModel3Run <- mxRun(myModel3)
328
329
330 First of all, let us view *myModel3* and more specifically the values of the matrices within that model.  Note that the ``$matrices`` lists one matrix, "Amatrix", and that the ``$algebras`` lists another, "Bmatrix".  To view values of matrices created with the ``mxMatrix()`` function, the argument is ``$values``; for matrices created with the ``mxAlgebra()`` function, the argument is ``$result``.  Note that when viewing a specific matrix, you can omit the ``$matrices`` or the ``$algebras`` arguments.
331
332    .. code-block:: r
333
334         >     myModel3
335         MxModel 'model3' 
336         type : default 
337         $matrices : 'Amatrix' 
338         $algebras : 'Bmatrix' 
339         $constraints :  
340         $intervals :  
341         $latentVars : none
342         $manifestVars : none
343         $data : NULL
344         $submodels :  
345         $expectation : NULL 
346         $fitfunction : NULL 
347         $compute : NULL 
348         $independent : FALSE 
349         $options :  
350         $output : FALSE 
351
352    .. code-block:: r
353
354         > myModel3$Amatrix$values
355              [,1]
356         [1,]    4
357
358    .. code-block:: r
359
360         > myModel3$Bmatrix$result
361         <0 x 0 matrix>
362
363 Fitted Model
364 ^^^^^^^^^^^^
365
366 Given we're looking at the model *myModel3* before it is run, results of algebra have not been computed yet.  Let us see how things change after running the model and viewing *myModel3Run*.
367
368    .. code-block:: r
369
370         >     myModel3Run
371         MxModel 'model3' 
372         type : default 
373         $matrices : 'Amatrix' 
374         $algebras : 'Bmatrix' 
375         $constraints :  
376         $intervals :  
377         $latentVars : none
378         $manifestVars : none
379         $data : NULL
380         $submodels :  
381         $expectation : NULL 
382         $fitfunction : NULL 
383         $compute : NULL 
384         $independent : FALSE 
385         $options :  
386         $output : TRUE
387
388    .. code-block:: r
389
390         > myModel3Run$Amatrix$values
391              [,1]
392         [1,]    4
393    
394    .. code-block:: r
395
396         > myModel3Run$Bmatrix$result
397              [,1]
398         [1,]    5
399
400 You will notice that the structure of the MxModel objects is identical, the value of the "Amatrix" has not changed, as it was a fixed element.  However, the value of the "Bmatrix" is now the result of the operation on the "Amatrix".  Note that we're here looking at the "Bmatrix" within the MxModel object *myModel3Run*.   Please verify that the original MxAlgebra objects *myBmatrix* and *myAmatrix* remain unchanged.  The mxModel() function call has made its own internal copies of these objects, and it is only these internal copies that are being manipulated.  In computer science terms, this is referred to as *pass by value*.
401
402
403 Pass By Value
404 ^^^^^^^^^^^^^
405
406 Let us insert a mini-lecture on the R programming language.  Our experience has found that this exercise will greatly increase your understanding of the OpenMx language. 
407
408 As this is such a crucial concept in R (unlike many other programming languages), let us look at it in a simple R example.  We will start by assigning the value 4 to the object "avariable", and then display it.  If we then add 1 to this object, and display it again, notice that the value of "avariable" has not changed.
409
410 .. cssclass:: input
411 ..
412
413    R Code
414    
415    .. code-block:: r
416
417         avariable <- 4
418         avariable
419     
420         avariable +1
421         avariable
422     
423 Now we introduce a function, as OpenMx is a collection of purposely built functions.  The function takes a single argument (the object "number"), adds one to the argument "number" and assigns the result to "number", and then returns the incremented number back to the user.  This function is given the name ``addone()``.  We then apply the function to the object "avariable", as well as display "avariable".  Thus, the objects "addone" and "avariable" are defined. The object assigned to "addone" is a function, while the value assigned to "avariable" is the number 4. 
424
425 .. cssclass:: input
426 ..
427
428    R Code
429    
430    .. code-block:: r
431
432         addone <- function(number) {
433             number <- number + 1
434             return(number)
435             }
436
437         addone(avariable)
438         avariable
439
440 Note that it may be prudent to use the ``print()`` function to display the results back to the user.  When R is run from a script rather than interactively, results will not be displayed unless the function ``print()`` is used as shown below.
441
442 .. cssclass:: input
443 ..
444
445    R Code
446    
447    .. code-block:: r
448
449         print(addone(avariable))
450         print(avariable)
451
452 What is the result of executing this code? Try it. The correct results are 5 and 4.  But why is the object "avariable" still 4, even after the ``addone()`` function was called? The answer to this question is that R uses pass-by-value function call semantics.
453
454 In order to understand pass-by-value semantics, we must understand the difference between *objects* and *values*. The *objects* declared in this example are "addone", "avariable", and "number".  The *values* refer to the things that are stored by the *objects*.  In programming languages that use pass-by-value semantics, at the beginning of a function call it is the *values* of the argument list that are passed to the function.  
455
456 The object "avariable" cannot be modified by the function ``addone()``.  If I wanted to update the value stored in the object, I would have needed to replace the expression as follows:
457
458 .. cssclass:: input
459 ..
460
461    R Code
462    
463    .. code-block:: r
464
465         print(avariable <- addone(avariable))
466         print(avariable)
467     
468 Try it.  The updated example prints out 5 and 5.  The lesson from this exercise is that the only way to update a object in a function call is to capture the result of the function call [#f1]_.  This lesson is sooo important that we'll repeat it:
469
470 *the only way to update an object in a function call is to capture the result of the function call.*
471
472 R has several built-in types of values that you are familiar with: numerics, integers, booleans, characters, lists, vectors, and matrices. In addition, R supports S4 object values to facilitate object-oriented programming.  Most of the functions in the OpenMx library return S4 object values.  You must always remember that R does not discriminate between built-in types and S4 object types in its call semantics.  Both built-in types and S4 object types are passed by value in R (unlike many other languages).
473
474 .. rubric:: Footnotes
475
476 .. [#f1] There are a few exceptions to this rule, but you can be assured such trickery is not used in the OpenMx library.
477
478
479 Styles
480 ------
481
482 In the beginning of the introduction, we discussed three styles of writing OpenMx code: the piecewise, stepwise and classic styles.  Let's take the most recent model and show how it can be written in these three styles.  
483
484 Piecewise Style
485 ^^^^^^^^^^^^^^^
486
487 The style we used in *myModel3* is the piecewise style.  We repeat the different statements here for clarity
488
489 .. cssclass:: input
490 ..
491
492    OpenMx Code
493    
494    .. code-block:: r
495
496         myAmatrix <- mxMatrix(type="Full", nrow=1, ncol=1, values=4, name="Amatrix")
497         myBmatrix <- mxAlgebra(expression=Amatrix+1, name="Bmatrix")
498
499         myModel3 <- mxModel(myAmatrix,myBmatrix, name="model3")
500     
501         myModel3Run <- mxRun(myModel3)
502
503 Each argument of the ``mxModel()`` statement is defined separately first as independent R objects which are then combined in one model statement.
504
505 Stepwise Style
506 ^^^^^^^^^^^^^^^
507
508 For the stepwise style, we start with an ``mxModel()`` with just one argument, as we originally did with the "Amatrix" in *myModel1*, as repeated below.  We could run this model to make sure it's syntactically correct.
509
510 .. cssclass:: input
511 ..
512
513     OpenMx Code
514
515     .. code-block:: r
516
517         myModel1 <- mxModel(
518                         mxMatrix(type="Full", nrow=1, ncol=1, values=4, name="Amatrix") 
519                     )
520
521         myModel1Run <- mxRun(myModel1)
522  
523 Then we would build a new model starting from the first model.  To do this, we invoke a special feature of the first argument of an ``mxModel()``.  If it is the name of a saved MxModel object, for example *myModel1*, the arguments of that model would be automatically included in the new model.  These arguments can be changed (or not) and new arguments can be added.  Thus, in our example, where we want to keep the "Amatrix" and add the "Bmatrix", our second model would look like this.  
524
525 .. cssclass:: input
526 ..
527
528    OpenMx Code
529    
530    .. code-block:: r
531
532         myModel4 <- mxModel(myModel1,
533                         mxAlgebra(expression=Amatrix+1, name="Bmatrix"), 
534                         name="model4"
535                         )
536         
537         myModel4Run <- mxRun(myModel4)
538     
539 Note that we call it "model4", by adding a ``name`` argument to the ``mxModel()`` as to not overwrite our previous "model1".
540
541 Classic Style
542 ^^^^^^^^^^^^^
543
544 The final style may be reminiscent of classic Mx.  Here we build all the arguments explicitly within one ``mxModel()``.  As a result only one R object is created prior to ``mxRun()`` ning the model.  This style is more compact than the others but harder to debug.
545
546 .. cssclass:: input
547 ..
548
549    OpenMx Code
550    
551    .. code-block:: r
552
553         myModel5 <- mxModel(
554                         mxMatrix(type="Full", nrow=1, ncol=1, values=4, name="Amatrix"),
555                         mxAlgebra(expression=Amatrix+1, name="Bmatrix"), 
556                         name="model5"
557                         )
558         
559         myModel5Run <- mxRun(myModel5)
560
561 You may have seen an alternative version with the first argument in quotes.  In that case, that argument refers to the name of the model and not to a previously defined model.  Thus, the following specification is identical to the previous one.  Note also that it is not necessary to add the 'names' of the arguments, as long as the arguments are listed in their default order, which can easily be verified by using the standard way to get help about a function (in this case ``?mxMatrix()`` ).
562
563 .. cssclass:: input
564 ..
565
566    OpenMx Code
567    
568    .. code-block:: r
569
570         myModel5 <- mxModel("model5",
571                         mxMatrix(type="Full", nrow=1, ncol=1, values=4, name="Amatrix"),
572                         mxAlgebra(expression=Amatrix+1, name="Bmatrix")
573                         )
574         
575         myModel5run <- mxRun(myModel5)
576
577 Note that all arguments are separated by commas.  In this case, we've also separated the arguments on different lines, but that is only for clarity.  No comma is needed after the last argument!  If you accidentally put one in, you get the generic error message *'argument is missing, with no default'* meaning that you forgot something and R doesn't know what it should be. The bracket on the following line closes the ``mxModel()`` statement.
578
579
580 Data functions
581 --------------
582
583 Most models will be fitted to data, not just a single number.  We will briefly introduce how to read data that are pre-packaged with the OpenMx library as well as reading in your own data.  All standard R utilities can be used here.  The critical part is to run an OpenMx model on these data, thus another OpenMx function ``mxData()`` is needed.
584
585 Reading Data
586 ^^^^^^^^^^^^
587
588 The ``data`` function can be used to read sample data that has been pre-packaged into the OpenMx library. One such sample data set is called "demoOneFactor".  
589
590 .. cssclass:: input
591 ..
592
593    R Code
594    
595    .. code-block:: r
596
597         data(demoOneFactor)
598
599 In order to read your own data, you will most likely use the ``read.table``, ``read.csv``, ``read.delim`` functions, or other specialized functions available from CRAN to read from 3rd party sources.  We recommend you install the package **psych** which provides succinct descriptive statistics with the ``describe()`` function.
600
601 .. cssclass:: input
602 ..
603
604    R Code
605    
606    .. code-block:: r
607    
608         require(psych)
609         describe(demoOneFactor)
610
611 The output of this function is shown below.
612
613     .. code-block:: r
614
615            var   n  mean   sd median trimmed  mad   min  max range  skew kurtosis   se
616         x1   1 500 -0.04 0.45  -0.03   -0.04 0.46 -1.54 1.22  2.77 -0.05     0.01 0.02
617         x2   2 500 -0.05 0.54  -0.03   -0.04 0.55 -2.17 1.72  3.89 -0.14     0.05 0.02
618         x3   3 500 -0.06 0.61  -0.03   -0.05 0.58 -2.29 1.83  4.12 -0.17     0.23 0.03
619         x4   4 500 -0.06 0.73  -0.08   -0.05 0.75 -2.48 2.45  4.93 -0.08     0.05 0.03
620         x5   5 500 -0.08 0.82  -0.08   -0.07 0.89 -2.62 2.18  4.80 -0.10    -0.23 0.04
621     
622
623 Now that the data are accessible in R, we need to make them readable into our OpenMx model.
624
625 Data Source 
626 ^^^^^^^^^^^
627
628 A ``mxData()`` function is used to construct a data source for the model.  OpenMx can handle fitting models to summary statistics and to raw data.
629
630 The most commonly used **summary statistics** are covariance matrices, means and correlation matrices; information on the variances is lost/unavailable with correlation matrices, so these are usually not recommended.
631
632 These days, the standard approach for model fitting applications is to use **raw data**, which is simply a data table or rectangular file with columns representing variables and rows representing subjects.  The primary benefit of this approach is that it handles datasets with missing values very conveniently and appropriately.
633
634 Covariance Matrix
635 ^^^^^^^^^^^^^^^^^
636
637 We will start with an example using summary data, so we are specifying a covariance matrix by using the R function ``cov`` to generate a covariance matrix from the data frame.  In addition to reading in the actual covariance matrix as the first (``observed``) argument, we specify the ``type`` (one of "cov","cor","sscp" and "raw") and the number of observations (``numObs``).
638
639 .. cssclass:: input
640 ..
641
642    OpenMx Code
643    
644    .. code-block:: r
645
646         exampleDataCov <- mxData(observed=cov(demoOneFactor), type="cov", numObs=500)
647     
648 We can view what *exampleDataCov* looks like for OpenMx.
649
650     .. code-block:: r
651
652          >      exampleDataCov
653          MxData 'data' 
654          type : 'cov' 
655          numObs : '500' 
656          Data Frame or Matrix : 
657                    x1        x2        x3        x4        x5
658          x1 0.1985443 0.1999953 0.2311884 0.2783865 0.3155943
659          x2 0.1999953 0.2916950 0.2924566 0.3515298 0.4019234
660          x3 0.2311884 0.2924566 0.3740354 0.4061291 0.4573587
661          x4 0.2783865 0.3515298 0.4061291 0.5332788 0.5610769
662          x5 0.3155943 0.4019234 0.4573587 0.5610769 0.6703023
663          Means : NA 
664          Acov : NA 
665          Thresholds : NA
666     
667 Some models may include predictions for the mean(s).  We could add an additional ``means`` argument to the ``mxData`` statement to read in the means as well.
668
669 .. cssclass:: input
670 ..
671
672    OpenMx Code
673    
674    .. code-block:: r
675
676         exampleDataCovMeans <- mxData(observed=cov(demoOneFactor), 
677                                    means=(colMeans(demoOneFactor), type="cov", numObs=500)
678     
679 The output for *exampleDataCovMeans* would have the following extra lines.
680
681     .. code-block:: r
682
683         ....
684         Means : 
685                       x1          x2          x3          x4          x5
686         [1,] -0.04007841 -0.04583873 -0.05588236 -0.05581416 -0.07555022
687         Acov : NA 
688         Thresholds : NA
689     
690 Raw Data
691 ^^^^^^^^
692
693 Note that for most real life examples, raw data are the preferred option, except in cases where complete data are available on all variables included in the analyses.  In that situation, using summary statistics is faster.  To change the current example to use raw data, we would read in the data explicitly and specify the ``type`` as "raw".  The ``numObs`` is no longer required as the sample size is counted automatically.
694
695 .. cssclass:: input
696 ..
697
698    OpenMx Code
699    
700    .. code-block:: r
701
702         exampleDataRaw <- mxData(observed=demoOneFactor, type="raw")
703
704 Printing this MxData object would result in listing the whole data set.  We show just the first few lines here:
705
706     .. code-block:: r
707
708          > exampleData
709          MxData 'data' 
710          type : 'raw' 
711          numObs : '500' 
712          Data Frame or Matrix : 
713                         x1            x2           x3           x4           x5
714          1   -1.086832e-01 -0.4669377298 -0.177839881 -0.080931127 -0.070650263
715          2   -1.464765e-01 -0.2782619339 -0.273882553 -0.154120074  0.092717293
716          3   -6.399140e-01 -0.9295294042 -1.407963429 -1.588974090 -1.993461644
717          4    2.150340e-02 -0.2552252972  0.097330513 -0.117444884 -0.380906486
718          5    ....
719
720
721 The data to be used for our example are now ready in either **covariance matrix** or **raw data** format.
722
723 Model functions
724 ---------------
725
726 We introduce here several new features by building a basic factor model to real data.  A useful tool to represent such a model is drawing a path diagram which is mathematically equivalent to equations describing the model.  If you're not familiar with the method of path analysis, we suggest you read one of the key reference books [LI1986]_.
727
728 ..[LI1986]  Li, C.C. (1986). Path Analysis - A Primer.  The Boxwood Press, Pacific Grove, CA.
729
730 Briefly, squares are used for observed variables; latent variables are drawn in circles.  One-headed arrows are drawn to represent causal relationships.  Correlations between variables are represented with two-headed arrows.  Double-headed paths are also used for variances of variables.  Below is a figure of a one factor model with five indicators (x1..x5). We have added a value of 1.0 to the variance of the latent variable **G** as a fixed value.  All the other paths in the models are considered free parameters and are to be estimated.
731
732 .. image:: graph/OneFactorModel.png
733     :height: 2in
734     
735 Variables
736 ^^^^^^^^^
737
738 To specify this path diagram in OpenMx, we need to indicate which variables are observed or manifest and which are latent.  The ``mxModel()`` arguments ``manifestVars`` and ``latentVars`` both take a vector of variable names.   In this case the manifest variables are "x1", "x2", "x3", "x4", "x5" and the latent variable is "G".  The R function ``c()`` is used to build the vectors.
739
740 .. cssclass:: input
741 ..
742
743    OpenMx Code
744    
745    .. code-block:: r
746
747         manifests <- c("x1","x2","x3","x4","x5")
748         latents <- c("G")
749         
750         manifestVars = manifests
751         latentVars = latents
752
753 This could be written more succinctly as follows.
754
755 .. cssclass:: input
756 ..
757
758    OpenMx Code
759    
760    .. code-block:: r
761    
762         manifestVars = names(demoOneFactor)
763         latentVars = c("G")
764
765 because the R ``names()`` function call returns the vector of names that we want (the observed variables in the data frame "demoOneFactor").
766
767 Path Creation
768 ^^^^^^^^^^^^^
769
770 Paths are created using the ``mxPath()`` function. Multiple paths can be created with a single invocation of the ``mxPath()`` function. 
771
772 - The ``from`` argument specifies the path sources, and the ``to`` argument specifies the path sinks.  If the ``to`` argument is missing, then it is assumed to be identical to the ``from`` argument. 
773 - The ``connect`` argument specifies the type of the source to sink connection, which can be one of five types.  For our example, we use the default "single" type in which the :math:`i^{th}` element of the ``from`` argument is matched with the :math:`i^{th}` element of the ``to`` argument, in order to create a path.  
774 - The ``arrows`` argument specifies whether the path is unidirectional (single-headed arrow, "1") or bidirectional (double-headed arrow, "2").  
775 - The next three arguments are vectors: ``free``, is a boolean vector that specifies whether a path is free or fixed; ``values`` is a numeric vector that specifies the starting value of the path; ``labels`` is a character vector that assigns a label to each free or fixed parameter.  Paths with the same labels are constrained to be equal, and OpenMx insists that paths equated in this way have the same fixed or free status; if this is not the case it will report an error.
776
777 To specify the path model above, we need to specify three different sets of paths.  The first are the single-headed arrows from the latent to the manifest variables, which we will put into the R object *causalPaths* as they represent causal paths.  The second set are the residuals on the manifest variables, referred to as *residualVars*.  The third ``mxPath()`` statement fixes the variance of the latent variable to one, and is called *factorVars*.
778
779 .. cssclass:: input
780 ..
781
782    OpenMx Code
783    
784    .. code-block:: r
785
786         causalPaths  <- mxPath(from=latents, to=manifests)
787         residualVars <- mxPath(from=manifests, arrows=2)
788         factorVars   <- mxPath(from=latents, arrows=2, free=FALSE, values=1.0)
789
790 Note that several arguments are optional.  For example, we omitted the ``free`` argument for *causalPaths* and *residualVars* because the default is 'TRUE' which applies in our example.  We also omitted the ``connect`` argument for all three paths.  The default "single" type automatically generates paths from every variable back to itself for all the variances, both the *residualVars* or the *factorVars*, as neither of those statements includes the ``to`` argument.  For the *causalPaths*, the default ``connect`` type will generate separate paths from the latent to each of the manifest variables.  To keep things simple, we did not include ``values`` or ``labels`` arguments as they are not strictly needed for this example, but this may not be true in general.  Once the variables and paths have been specified, the predicted covariance matrix will be generated from the implied path diagram in the backend of OpenMx using the RAM notation (see below).
791
792 Equations
793 ^^^^^^^^^
794
795 For those more in tune with equations and matrix algebra, we can represent the model using matrix algebra rather than path specifications.  For reasons that may become clear later, the expression for the expected covariances between the manifest variables is given by  
796
797 .. math::
798    :nowrap:
799
800    \begin{eqnarray*} 
801    \mbox{Cov} ( x_{ij}) = facLoadings * facVariances * facLoadings^\prime + resVariances
802    \end{eqnarray*}
803
804 where *facLoadings* is a column vector of factor loadings, *facVariances* is a symmetric matrix of factor variances and *resVariances* is a diagonal matrix of residual variances.  You might have noticed the correspondence between *causalPaths* and *facLoadings*, between *residualVars* and *resVariances*, and between *factorVars* and *facVariances*.  To translate this model into OpenMx using the matrix specification, we will define the three matrices first using the ``mxMatrix()`` function, and then specify the algebra using the ``mxAlgebra()`` function.
805
806 Matrix Creation
807 ^^^^^^^^^^^^^^^
808
809 The next three lines create three ``MxMatrix()`` objects, using the ``mxMatrix()`` function.  The first argument declares the ``type`` of the matrix, the second argument declares the number of rows in the matrix (``nrow``), and the third argument declares the number of columns (``ncol``).  The ``free`` argument specifies whether an element is a free or fixed parameter.  The ``values`` argument specifies the starting values for the elements in the matrix, and the ``name`` argument specifies the name of the matrix. 
810
811 .. cssclass:: input
812 ..
813
814    OpenMx Code
815    
816    .. code-block:: r
817
818         mxFacLoadings  <-  mxMatrix(type="Full", nrow=5, ncol=1, 
819                                     free=TRUE, values=0.2, name="facLoadings")
820         mxFacVariances <-  mxMatrix(type="Symm", nrow=1, ncol=1, 
821                                     free=FALSE, values=1, name="facVariances")
822         mxResVariances <-  mxMatrix(type="Diag", nrow=5, ncol=5, 
823                                     free=TRUE, values=1, name="resVariances")
824
825 Each ``MxMatrix()`` object is a container that stores five matrices of equal dimensions.  The five matrices stored in a ``MxMatrix()`` object are:``free``, ``values``, ``labels``, ``lbound``, and ``ubound``.  ``Free`` stores a boolean vector that determines whether a element is free or fixed.  ``Values`` stores the current values of each element in the matrix.  ``Labels`` stores a character label for each element in the matrix. And ``lbound`` and ``ubound`` store the lower and upper bounds, respectively, for each element that is a free parameter.  If a element has no label, lower bound, or upper bound, then an NA value is stored in the element of the respective matrix.
826  
827 Algebra Creation
828 ^^^^^^^^^^^^^^^^
829
830 An ``mxAlgebra()`` function is used to construct an expression for any algebra, in this case the expected covariance algebra.  The first argument (``expression``) is the algebra expression that will be evaluated by the numerical optimizer.  The matrix operations and functions that are permitted in an MxAlgebra expression are listed in the help for the ``mxAlgebra`` function (obtained by ``?mxAlgebra``).  The algebra expression refers to entities according to ``name`` argument of the MxMatrix objects.
831
832 .. cssclass:: input
833 ..
834
835    OpenMx Code
836    
837    .. code-block:: r
838
839         mxExpCov <- mxAlgebra(expression=facLoadings %*% facVariances %*% t(facLoadings) 
840                                         + resVariances, name="expCov")
841
842 You can see a direct correspondence between the formula above and the expression used to create the expected covariance matrix *myExpCov*.
843
844 Expectation-FitFunction
845 -----------------------
846
847 To fit a model to data, the differences between the observed covariance matrix (the data, in this case the summary statistics) and model-implied expected covariance matrix are minimized using a fit function.  Fit functions are functions for which free parameter values are chosen such that the value of the fit function is minimized.  Now that we have specified data objects and path or matrix/algebra objects for the predicted covariances of our model, we need to link the two and execute them which is typically done with ``Expectation`` and ``FitFunction`` statements.  PS. These two statements replace the mxObjective functions in earlier versions of OpenMx.  
848
849 RAM Expectation 
850 ^^^^^^^^^^^^^^^
851
852 When using a path specification of the model, the fit function is always ``RAM`` which is indicated by using the ``type`` argument.  We don't have to specify the fit function explicitly with an ``mxExpectation()`` and ``FitFunction`` argument, instead we simply add the following argument to the model.
853
854 .. cssclass:: input
855 ..
856
857    OpenMx Code
858    
859    .. code-block:: r
860
861         type="RAM"
862     
863 To gain a better understanding of the RAM principles, we recommend reading [RAM1990]
864
865 ..[RAM1990]  McArdle, J.J. & Boker, S.M. (1990). RAMpath: Path diagram software. Denver: Data Transforms Inc.
866
867 Normal Expectation
868 ^^^^^^^^^^^^^^^^^^
869
870 When using a matrix specification, ``mxExpectationNormal()`` defines how model expectations are calculated using the matrices/algebra implied by the ``covariance`` argument and optionally the ``means``.  For this example, we are specifying an expected covariance algebra (``covariance``) omitting an expected means algebra.  The expected covariance algebra is referenced according to its name, i.e. the ``name`` argument of the MxAlgebra created above.  We also need to assign ``dimnames`` for the rows and columns of this covariance matrix, such that a correspondence can be determined between the expected and the observed covariance matrices.  Subsequently we are specifying a maximum likelihood fit function with the ``mxFitFunctionML()`` statement.
871
872 .. cssclass:: input
873 ..
874
875    OpenMx Code
876    
877    .. code-block:: r
878
879         expectCov <- mxExpectationNormal(covariance="expCov", 
880                                          dimnames=names(demoOneFactor))
881         funML <- mxFitFunctionML()
882
883 The above expectation and fit function can be used when fitting to covariance matrices.  A model for the predicted means is optional.  However, when fitting to raw data, an expectation has to be used that specifies both a model for the means and for the covariance matrices, paired with the appropriate fit function.  In the case of raw data, the ``mxFitFunctionML()`` function uses full-information maximum likelihood to provide maximum likelihood estimates of free parameters in the algebra defined by the ``covariance`` and ``means`` arguments. The ``covariance`` argument takes an ``MxMatrix`` or ``MxAlgebra`` object, which defines the expected covariance of an associated ``MxData`` object. Similarly, the ``means`` argument takes an ``MxMatrix`` or ``MxAlgebra`` object to define the expected means of an associated ``MxData`` object. The ``dimnames`` arguments takes an optional character vector. This vector is assigned to be the ``dimnames`` of the means vector, and the row and columns ``dimnames`` of the covariance matrix. 
884
885 .. cssclass:: input
886 ..
887
888    OpenMx Code
889    
890    .. code-block:: r
891
892         expectCovMeans <- mxExpectationNormal(covariance="expCov", means="expMeans", 
893                                               dimnames=names(demoOneFactor))
894         funML <- mxFitFunctionML()
895
896 Raw data can come in two forms, continuous or categorical.  While **continuous data** have an unlimited number of possible values, their frequencies typically form a normal distribution.
897
898 There are basically two flavors of **categorical data**.  If only two response categories exist, for example Yes and No, or affected and unaffected, we are dealing with binary data.  Variables with three or more categories are considered ordinal.
899
900 Continuous Data
901 ^^^^^^^^^^^^^^^
902
903 When the data to be analyzed are continuous, and models are fitted to raw data, the ``mxFitFunctionML()`` function will take two arguments, the ``covariance`` and the ``means`` argument, as well as ``dimnames`` to match them up with the observed data.
904
905 .. cssclass:: input
906 ..
907
908    OpenMx Code
909    
910    .. code-block:: r
911
912         expectRaw <- mxExpectationNormal(covariance="expCov", means="expMeans", 
913                                          dimnames=manifests)
914         funML <- mxFitFunctionFIML()
915
916 If the variables to be analyzed have at least 15 possible values, we recommend to treat them as continuous data.  As will be discussed later in the documentation, the power of the study is typically higher when dealing with continuous rather than categorical data.
917
918 Categorical Data
919 ^^^^^^^^^^^^^^^^
920
921 For categorical - be they binary or ordinal - data, an additional argument is needed for the ``mxFitFunctionML()`` function, besides the ``covariance`` and ``means`` arguments, namely the ``thresholds`` argument.
922     
923 .. cssclass:: input
924 ..
925
926     OpenMx Code
927
928     .. code-block:: r
929
930         expFunOrd <- mxExpectationNormal(covariance="expCov", means="expMeans", 
931                                          thresholds="expThres", dimnames=manifests)
932         funML <- mxFitFunctionFIML()
933
934 For now, we will stick with the factor model example and fit it to covariance matrices, calculated from the raw continuous data.
935
936
937 Methods
938 -------
939
940 We have introduced two ways to create a model.  One is the **path method**, in which observed and latent variables are specified as well as the causal and correlational paths that connect the variables to form the model.  This method may be more intuitive as the model maps on directly to the diagram.  This of course assumes that the path diagram is drawn mathematically correct.  Once the model is 'drawn' or specified correctly in this way, OpenMx translates the paths into RAM notation for predicted covariance matrices.
941
942 Alternatively, we can specify the model using the **matrix method** by creating the necessary matrices and combining them using algebra to generate the expected covariance matrices (and optionally the mean/threshold vectors).  Although less intuitive, this method provides greater flexibility for developing more complex models.  Let us look at examples of both.
943
944 Path Method
945 ^^^^^^^^^^^
946
947 We have previously generated all the pieces that go into the model, using the path method specification.  As we have discussed before, the ``mxModel()`` function is somewhat of a swiss-army knife.  The first argument to the ``mxModel()`` function can be an argument of type ``name`` (and appear in quotes), in which case it is a newly generated model, or it can be a previously defined model object.  In the latter case, the new model 'inherits' all the characteristics (arguments) of the old model, which can be changed with additional arguments.  An ``mxModel()`` can contain ``mxData()``, ``mxPath()``, ``mxExpectation()``, ``mxFitFunction`` and other ``mxModel()`` statements as arguments.
948
949 The  following ``mxModel()`` function is used to create the 'one-factor' model, shown on the path diagram above.  The first argument is a ``name``, thus we are specifying a new model, called "One Factor".  By specifying the ``type`` argument to equal "RAM", we create a path style model. A RAM style model must include a vector of manifest variables (``manifestVars=``) and a vector of latent variables (``latentVars=``).   We then include the arguments for reading the example data *exampleDataCov*, and those that specify the paths of the path model *causalPaths*, *residualVars*, and *factorVars* which we created previously.
950
951 .. cssclass:: input
952 ..
953
954    OpenMx Code
955    
956    .. code-block:: r
957
958         factorModel1 <- mxModel(name="One Factor", 
959             type="RAM",
960             manifestVars = manifests,
961             latentVars = latents,
962             exampleDataCov, causalPaths, residualVars, factorVars)
963
964 When we display the contents of this model, note that we now have manifest and latent variables specified.  By using ``type`` ="RAM" we automatically use the expectation ``mxExpectationRAM`` which translates the path model into RAM specification [RAM1990] as reflected in the matrices **A**, **S** and **F**,  and the function ``mxFitFunctionML()``.  Briefly, the **A** matrix contains the asymmetric paths, which are the unidirectional paths in the *causalPaths* object, and represent the factor loadings from the latent variable onto the manifest variables.  The **S** matrix contains the symmetric paths which include both the bidirectional paths in *residualVars* and in *factorVars*.  The **F** matrix is the filter matrix.
965
966 The formula :math:`F(I-A)~*S*(I-A)~'F'`, where I is an identity matrix, ~ denotes the inverse and ' the transpose, generates the expected covariance matrix.
967
968    .. code-block:: r
969
970         >     factorModel1
971         MxModel 'One Factor' 
972         type : RAM 
973         $matrices : 'A', 'S', and 'F' 
974         $algebras :  
975         $constraints :  
976         $intervals :  
977         $latentVars : 'G' 
978         $manifestVars : 'x1', 'x2', 'x3', 'x4', and 'x5' 
979         $data : 5 x 5 
980         $data means : NA
981         $data type: 'cov' 
982         $submodels :  
983         $expectation : MxExpectationRAM 
984         $fitfunction : MxFitFunctionML 
985         $compute : NULL 
986         $independent : FALSE 
987         $options :  
988         $output : FALSE 
989
990 You can verify that after running the model, the new R object *factorFit* has similar arguments, except that they now contain the estimates from the model rather than the starting values.  For example, we can look at the values in the **A** matrix in the built model *factorModel*, and in the fitted model *factorFit*.  We will get back to this later.  Note also that from here on out, we use the convention the R object containing the built model will end with *Model* while the R object containing the fitted model will end with *Fit*.
991
992 .. cssclass:: input
993 ..
994
995    OpenMx Code
996    
997    .. code-block:: r
998
999         factorFit1 <- mxRun(factorModel1)
1000
1001 We can inspect the values of the **A** matrix in *factorModel1* and *factorFit1* respectively as follows.
1002
1003     .. code-block:: r
1004
1005         > factorModel1$A$values
1006            x1 x2 x3 x4 x5 G
1007         x1  0  0  0  0  0 0
1008         x2  0  0  0  0  0 0
1009         x3  0  0  0  0  0 0
1010         x4  0  0  0  0  0 0
1011         x5  0  0  0  0  0 0
1012         G   0  0  0  0  0 0
1013
1014         > factorFit1$A$values 
1015            x1 x2 x3 x4 x5         G
1016         x1  0  0  0  0  0 0.3971521
1017         x2  0  0  0  0  0 0.5036611
1018         x3  0  0  0  0  0 0.5772414
1019         x4  0  0  0  0  0 0.7027737
1020         x5  0  0  0  0  0 0.7962500
1021         G   0  0  0  0  0 0.0000000
1022
1023 We can also specify all the arguments directly within the ``mxModel()`` function, using the **classical** style, as follows.  The script reads data from disk, creates the one factor model, fits the model to the observed covariances, and prints a summary of the results. 
1024
1025 .. cssclass:: input
1026 ..
1027
1028    OpenMx Code
1029    
1030    .. code-block:: r
1031
1032         data(demoOneFactor)
1033         manifests <- names(demoOneFactor)
1034         latents <- c("G")
1035         
1036         factorModel1 <- mxModel(name="One Factor", 
1037             type="RAM",
1038             manifestVars=manifests,
1039             latentVars=latents,
1040             mxPath(from=latents, to=manifests),
1041             mxPath(from=manifests, arrows=2),
1042             mxPath(from=latents, arrows=2, free=FALSE, values=1.0), 
1043             mxData(observed=cov(demoOneFactor), type="cov", numObs=500)
1044         )
1045         
1046         factorFit1 <- mxRun(factorModel1)
1047         summary(factorFit1)
1048     
1049 For more details about the summary and alternative options to display model results, see below.
1050
1051 Matrix Method
1052 ^^^^^^^^^^^^^
1053
1054 We will now re-create the model from the previous section, but this time we will use a matrix specification technique. The script reads data from disk, creates the one factor model, fits the model to the observed covariances, and prints a summary of the results. 
1055
1056 We have already created separate objects for each of the parts of the model, which can then be combined in an ``mxModel()`` statement at the end.  To repeat ourselves, the name of an OpenMx entity bears no relation to the R object that is used to identify the entity. In our example, the object "mxFacLoadings" stores a value that is a MxMatrix object with the name "facLoadings".  Note, however, that it is not necessary to use different names for the name within the ``mxMatrix`` object and the name of the R object generated with the statement.  For more complicated models, using the same name for both rather different entities, may make it easier to keep track of the various pieces.  For now, we will use different names to highlight which one should be used in which context.
1057  
1058 .. cssclass:: input
1059 ..
1060
1061    OpenMx Code
1062    
1063    .. code-block:: r
1064
1065         data(demoOneFactor)
1066         
1067         factorModel2 <- mxModel(name="One Factor",
1068             exampleDataCov, mxFacLoadings, mxFacVariances, mxResVariances, 
1069             mxExpCov, expectCov, funML)
1070         factorFit2 <- mxRun(factorModel2)
1071         summary(factorFit2)
1072
1073
1074 Alternatively, we can write the script in the **classical** style and specify  all the matrices, algebras, objective function and data as arguments to the ``mxModel()``.
1075
1076 .. cssclass:: input
1077 ..
1078
1079    OpenMx Code
1080    
1081    .. code-block:: r
1082
1083         data(demoOneFactor)
1084         
1085         factorModel2 <- mxModel(name="One Factor",
1086             mxMatrix(type="Full", nrow=5, ncol=1, free=TRUE, values=0.2, name="facLoadings"),
1087             mxMatrix(type="Symm", nrow=1, ncol=1, free=FALSE, values=1, name="facVariances"),
1088             mxMatrix(type="Diag", nrow=5, ncol=5, free=TRUE, values=1, name="resVariances"),
1089             mxAlgebra(expression=facLoadings %*% facVariances %*% t(facLoadings) 
1090                                 + resVariances, name="expCov"),
1091             mxExpectationNormal(covariance="expCov", dimnames=names(demoOneFactor)),
1092             mxFitFunctionML()
1093             mxData(observed=cov(demoOneFactor), type="cov", numObs=500)
1094         )
1095         
1096         factorFit2 <- mxRun(factorModel2)
1097         summary(factorFit2)
1098
1099 Now that we've specified the model with both methods, we can run both examples and verify that they indeed provide the same answer by inspecting the two fitted R objects "factorFit1" and "factorFit2".
1100
1101
1102 Output
1103 ------
1104
1105 We can generate output in a variety of ways.  As you might expect, the **summary** function summarizes the model, including data, model parameters, goodness-of-fit and run statistics.
1106
1107 Note that the fitted model is an R object that can be further manipulated, for example, to output specific parts of the model or to use it as a basis for developing an alternative model.
1108
1109 Model Summary
1110 ^^^^^^^^^^^^^
1111
1112 The summary function (``summary(modelname)``) is a convenient method for displaying the highlights of a model after it has been executed.  Many R functions have an associated ``summary()`` function which summarizes all key aspects of the model.  In the case of OpenMx, the ``summary(model)`` includes a summary of the data, a list of all the free parameters with their name, matrix element locators, estimate and standard error, as well as lower and upper bounds if those were assigned.  Currently the list of goodness-of-fit statistics printed include the number of observed statistics, the number of estimated parameters, the degrees of freedom, minus twice the log-likelihood of the data, the number of observations, the chi-square and associated p-value and several information criteria.  Various time-stamps and the OpenMx version number are also displayed.
1113
1114    .. code-block:: r
1115
1116         >     summary(factorFit1)
1117         data:
1118         $`One Factor.data`
1119         $`One Factor.data`$cov
1120                   x1        x2        x3        x4        x5
1121         x1 0.1985443 0.1999953 0.2311884 0.2783865 0.3155943
1122         x2 0.1999953 0.2916950 0.2924566 0.3515298 0.4019234
1123         x3 0.2311884 0.2924566 0.3740354 0.4061291 0.4573587
1124         x4 0.2783865 0.3515298 0.4061291 0.5332788 0.5610769
1125         x5 0.3155943 0.4019234 0.4573587 0.5610769 0.6703023
1126
1127
1128         free parameters:
1129           name matrix row col   Estimate   Std.Error Std.Estimate      Std.SE lbound ubound
1130         1  One Factor.A[1,6]      A  x1   G 0.39715182 0.015549708   0.89130932 0.034897484              
1131         2  One Factor.A[2,6]      A  x2   G 0.50366066 0.018232433   0.93255458 0.033758321              
1132         3  One Factor.A[3,6]      A  x3   G 0.57724092 0.020448313   0.94384664 0.033435037              
1133         4  One Factor.A[4,6]      A  x4   G 0.70277323 0.024011318   0.96236250 0.032880581              
1134         5  One Factor.A[5,6]      A  x5   G 0.79624935 0.026669339   0.97255562 0.032574489              
1135         6  One Factor.S[1,1]      S  x1  x1 0.04081418 0.002812716   0.20556770 0.014166734              
1136         7  One Factor.S[2,2]      S  x2  x2 0.03801997 0.002805791   0.13034196 0.009618951              
1137         8  One Factor.S[3,3]      S  x3  x3 0.04082716 0.003152305   0.10915353 0.008427851              
1138         9  One Factor.S[4,4]      S  x4  x4 0.03938701 0.003408870   0.07385841 0.006392303              
1139         10 One Factor.S[5,5]      S  x5  x5 0.03628708 0.003678556   0.05413557 0.005487924              
1140
1141         observed statistics:  15 
1142         estimated parameters:  10 
1143         degrees of freedom:  5 
1144         -2 log likelihood:  -3648.281 
1145         saturated -2 log likelihood:  -3655.665 
1146         number of observations:  500 
1147         chi-square:  7.384002 
1148         p:  0.1936117 
1149         Information Criteria: 
1150              df Penalty Parameters Penalty Sample-Size Adjusted
1151         AIC:  -2.615998           27.38400                   NA
1152         BIC: -23.689038           69.53008             37.78947
1153         CFI: 0.9993583 
1154         TLI: 0.9987166 
1155         RMSEA:  0.03088043 
1156         timestamp: 2014-04-10 10:23:07 
1157         frontend time: 0.02934313 secs 
1158         backend time: 0.005492926 secs 
1159         independent submodels time: 1.907349e-05 secs 
1160         wall clock time: 0.03485513 secs 
1161         cpu time: 0.03485513 secs 
1162         openmx version number: 999.0.0-3297 
1163          
1164
1165 The table of free parameters requires a little more explanation.  First, ``<NA>`` are given for the name of elements that were not assigned a label.  Second, the columns 'row' and 'col' display the variables at the tail of the paths and the variables at the head of the paths respectively.  Third, standard errors are calculated.  We will discuss the use of standard errors versus confidence intervals later on.
1166
1167 Model Evaluation 
1168 ^^^^^^^^^^^^^^^^
1169
1170 The ``mxEval()`` function should be your primary tool for observing and manipulating the final values stored within a MxModel object.  The simplest form of the ``mxEval()`` function takes two arguments: an ``expression`` and a ``model``. The expression can be **any** arbitrary expresssion to be evaluated in R.  That expression is evaluated, but the catch is that any named entities or parameter names are replaced with their current values from the model.  The model can be either a built or a fitted model.
1171
1172 .. cssclass:: input
1173 ..
1174
1175    OpenMx Code
1176    
1177    .. code-block:: r
1178
1179         myModel6 <- mxModel('topmodel', 
1180             mxMatrix('Full', 1, 1, values=1, free=TRUE, labels='p1', name='A'),
1181             mxModel('submodel', 
1182                 mxMatrix('Full', 1, 1, values=2, free=FALSE, labels='p2', name='B')
1183             )
1184         )
1185
1186         myModel6Run <- mxRun(myModel6)
1187         mxEval(A + submodel.B + p1 + p2, myModel6)       # initial values
1188         mxEval(A + submodel.B + p1 + p2, myModel6Run)    # final values
1189
1190 Note that the ``expression`` can include both matrices, algebras as well as matrix element labels, each taking on the value of the model specified in the ``model`` argument.  To reinforce an earlier point, it is not necessary to restrict the expression only to valid MxAlgebra expressions.  In the following example, we use the ``harmonic.mean()`` function from the psych package.
1191
1192 .. cssclass:: input
1193 ..
1194
1195    OpenMx Code
1196    
1197    .. code-block:: r
1198
1199         library(psych)
1200         nVars <- 4
1201         mxEval(nVars * harmonic.mean(c(A, submodel.B)), myModel6)
1202
1203 When the name of an entity in a model collides with the name of a built-in or user-defined function in R, the named entity will supercede the function.  We strongly advice against naming entities with the same name as the predefined functions or values in R, such as `c`, `T`, and `F` among others.
1204
1205 The ``mxEval()`` function allows the user to inspect the values of named entities without explicitly poking at the internals of the components of a model.  We encourage the use of ``mxEval()`` to look at the state of a model either before the execution of a model or after model execution.
1206
1207
1208 Indexing Operator
1209 ^^^^^^^^^^^^^^^^^
1210
1211 MxModel objects support the ``$`` operator, also known as the list indexing operator, to access all the components contained within a model.  Here is an example collection of models that will help explain the uses of the ``$`` operator:
1212
1213 .. cssclass:: input
1214 ..
1215
1216    OpenMx Code
1217    
1218    .. code-block:: r
1219    
1220         myModel7 <- 
1221             mxModel('topmodel', 
1222                 mxMatrix(type='Full', nrow=1, ncol=1, name='A'),
1223                 mxAlgebra(A, name='B'),
1224                 mxModel('submodel1', 
1225                     mxConstraint(topmodel1.A == topmodel1.B, name = 'C'),
1226                     mxModel('undersub1', 
1227                         mxData(diag(3), type='cov', numObs=10)
1228                     )
1229                 ),
1230                     mxModel('submodel2', 
1231                         mxFitFunctionAlgebra('topmodel1.A')
1232                     )
1233             )
1234
1235 Access Elements
1236 ^^^^^^^^^^^^^^^
1237
1238 The first useful trick is entering the string ``model$`` in the R interpreter and then pressing the TAB key.  You should see a list of all the named entities contained within the ``model`` object.
1239
1240     .. code-block:: r
1241
1242         > model$
1243         model$A                    
1244         model$B                    
1245         model$submodel1
1246         model$submodel2            
1247         model$submodel1.C          
1248         model$undersub1
1249         model$undersub1.data
1250         model$submodel2.fitfunction
1251
1252 The named entities of the model are displayed in one of three modes. 
1253
1254 #. All of the submodels contained within the parent model are accessed by using their unique model name (``submodel1``, ``submodel2``, and ``undersub1``).  
1255
1256 #. All of the named entities contained within the parent model are displayed by their names (``A`` and ``B``).  
1257
1258 #. All of the named entities contained by the submodels are displayed in the ``modelname.entityname`` format (``submodel1.C``, ``submodel2.objective``, and ``undersub1.data``). 
1259
1260 Modify Elements
1261 ^^^^^^^^^^^^^^^
1262
1263 The list indexing operator can also be used to modify the components of an existing model. There are three modes of using the list indexing operator to perform modifications, and they correspond to the three modes for accessing elements.
1264
1265 In the first mode, a submodel can be replaced using the unique name of the submodel or even eliminated.
1266
1267 .. cssclass:: input
1268 ..
1269
1270    OpenMx Code
1271    
1272    .. code-block:: r
1273
1274         # replace 'submodel1' with the contents of the mxModel() expression
1275         model$submodel1 <- mxModel(...)      
1276         # eliminate 'undersub1' and all children models
1277         model$undersub1 <- NULL              
1278
1279 In the second mode, the named entities of the parent model are modified using their names.  Existing matrices can be eliminated or new matrices can be created.
1280
1281 .. cssclass:: input
1282 ..
1283
1284    OpenMx Code
1285    
1286    .. code-block:: r
1287    
1288         # eliminate matrix 'A'
1289         model$A <- NULL
1290         # create matrix 'D'
1291         model$D <- mxMatrix(...)             
1292
1293 In the third mode, named entities of a submodel can be modified using the ``modelname.entityname`` format.  Again existing elements can be eliminated or new elements can be created.
1294
1295 .. cssclass:: input
1296 ..
1297
1298    OpenMx Code
1299    
1300    .. code-block:: r
1301    
1302         # eliminate constraint 'C' from submodel1
1303         model$submodel1.C <- NULL
1304         # create algebra 'D' in undersub1
1305         model$undersub1.D <- mxAlgebra(...)         
1306         # create 'undersub2' as a child model of submodel1
1307         model$submodel1.undersub2 <- mxModel(...)   
1308
1309 Keep in mind that when using the list indexing operator to modify a named entity within a model, the name of the created or modified entity is always the name on the left-hand side of the ``<-`` operator.  This feature can be convenient, as it avoids the need to specify a name of the entity on the right-hand side of the ``<-`` operator.
1310
1311
1312 Classes
1313 -------
1314
1315 We have introduced a number of OpenMx functions which correspond to specific classes which are summarized below. 
1316 The basic unit of abstraction in the OpenMx library is the model.  A model serves as a container for a collection of matrices, algebras, constraints, expectation, fit functions, data sources, and nested sub-models.  In the parlance of R, a model is a value that belongs to the class MxModel that has been defined by the OpenMx library.  The following table indicates what classes are defined by the OpenMx library.
1317
1318                     +--------------------+---------------------+
1319                     | entity             | S4 class            |
1320                     +====================+=====================+
1321                     | model              | MxModel             | 
1322                     +--------------------+---------------------+
1323                     | data source        | MxData              |
1324                     +--------------------+---------------------+
1325                     | matrix             | MxMatrix            |
1326                     +--------------------+---------------------+
1327                     | algebra            | MxAlgebra           |
1328                     +--------------------+---------------------+
1329                     | expectation        | MxExpectationRAM    |
1330                     |                    | MxExpectationNormal |
1331                     +--------------------+---------------------+
1332                     | fit function       | MxFitFunctionML     |
1333                     +--------------------+---------------------+                    
1334                     | constraint         | MxConstraint        |
1335                     +--------------------+---------------------+
1336
1337 All of the entities listed in the table are identified by the OpenMx library by the name assigned to them.  A name is any character string that does not contain the "." character.  In the parlance of the OpenMx library, a model is a container of named entities.  The name of an OpenMx entity bears no relation to the R object that is used to identify the entity. In our example, the object ``factorModel`` is created with the ``mxModel()`` function and stores a value that is an "MxModel" object with the name 'One Factor'.