A Linear Structural Equation Model for Covert Verb Generation Based on Independent Component Analysis of fMRI Data from Children and Adolescents

Human language is a complex and protean cognitive ability. Young children, following well defined developmental patterns learn language rapidly and effortlessly producing full sentences by the age of 3 years. However, the language circuitry continues to undergo significant neuroplastic changes extending well into teenage years. Evidence suggests that the developing brain adheres to two rudimentary principles of functional organization: functional integration and functional specialization. At a neurobiological level, this distinction can be identified with progressive specialization or focalization reflecting consolidation and synaptic reinforcement of a network (Lenneberg, 1967; Muller et al., 1998; Berl et al., 2006). In this paper, we used group independent component analysis and linear structural equation modeling (McIntosh and Gonzalez-Lima, 1994; Karunanayaka et al., 2007) to tease out the developmental trajectories of the language circuitry based on fMRI data from 336 children ages 5–18 years performing a blocked, covert verb generation task. The results are analyzed and presented in the framework of theoretical models for neurocognitive brain development. This study highlights the advantages of combining both modular and connectionist approaches to cognitive functions; from a methodological perspective, it demonstrates the feasibility of combining data-driven and hypothesis driven techniques to investigate the developmental shifts in the semantic network.

Moreover, a priori model testing based on fewer free parameters presents significant challenges for an investigation that deals with the developmental trajectories of skills sub-serving verb generation. In particular, such limitations significantly impact our ability to test specific models (e.g., regionally weighted or focal network models) of language development using neuroimage data as input to connectionist approaches for neurocognitive modeling. In the connectionist approaches a system behavior is captured by adjusting the weights on connections between elements in the network to investigate how the statistical structure of inputs influences the behavior of the network (Plaut et al., 1996). Therefore, with more parameters (degrees of freedom) one is better positioned to capture any development shifts in neurocognitive modeling.
The utility of independent component analysis (ICA) for examining changes in brain networks associated with age and brain development has recently been demonstrated in the context of resting state (Stevens et al., 2009a) as well as active neurocognitive processes (Stevens et al., 2009a,b) such as language function Karunanayaka et al., 2007Karunanayaka et al., , 2010Karunanayaka et al., , 2011Kim et al., 2011). Unlike model-based approaches, ICA is a datadriven technique capable of detecting additional task-related neural networks that exhibit activity with different temporal behavior (Calhoun et al., 2001a). This approach has significant advantages IntroductIon Functional brain imaging methods have recently emerged as means of investigating connectivity and the dynamic flow of information across neural networks sub-serving cognitive functions (McIntosh and Gonzalez-Lima, 1994;Friston et al., 2003;Penny et al., 2004a,b). These methods measure, e.g., the generated electrical/magnetic fields (EEG/MEG) or the hemodynamic response associated with neural activity (fMRI). The functional data analysis methods frequently focus on identifying areas of activation under different behavioral conditions with less attention paid to the behavior of the underlying network .
Until recently, fMRI studies have employed model-based approaches predicated upon a priori knowledge of an applied stimulus and the brain's response [hemodynamic response function (HRF)] to the stimulus (Bandettini et al., 1993;Worsley and Friston, 1995). Such models are typically based on canonical forms for the HRF and do not reflect individual variations or account for differences between individuals of different age, sex, or pathologies. We have previously discussed that this statistical approach may not capture the complexity of brain networks supporting a language task such as covert verb generation . Several methods have been proposed to circumvent this drawback by avoiding assumptions about the shape of the HRF (Ollinger et al., 2001).
A linear structural equation model for covert verb generation based on independent component analysis of fMRI data from children and adolescents prenatal period (Wada et al., 1975;Chi et al., 1977), which suggests a structural basis for early left hemisphere lateralization of related functions (Foundas et al., 1994). While functional asymmetries are not present at birth (Kotilahti et al., 2010) the anatomical asymmetries, in association with genetic factors, may underlie later development of functional asymmetries (Szaflarski et al., 2002;Francks et al., 2007). In fact, previous reports from the parent project that generated this data set (Szaflarski et al., 2006a,b;Holland et al., 2007) and others (Wood et al., 2004;Chou et al., 2006) have indicated that the initial left lateralization in this area strengthens with age with maximum left lateralization achieved around the age of 20-25 years followed by gradual decrease in the observed asymmetries with increasing age (Szaflarski et al., 2002(Szaflarski et al., , 2006a. Therefore, in this study we expected to confirm the age-related changes in the networks (inter and intra) that sub-serve verbal abilities.
Because the noun must be held in working memory as the verbs are generated, we expect that the temporal cortex must be connected to a fronto-parietal network that is routinely activated in studies involving working memory (Chein et al., 2003). This includes activation of the inferior frontal gyrus, dorsolateral prefrontal cortex, and parietal cortex. More specifically, the superior temporal cortex is functionally connected to the inferior frontal gyrus (through the arcuate fasciculus) as the anatomical connections between these two regions are well established (Catani et al., 2005). The dorsolateral prefrontal (executive control) and parietal (sustained attention to words) cortices modulate activity in this region through either the superior branch of the arcuate fasciculus (Catani et al., 2005) or the superior longitudinal fasciculus. Because working memory shows age-related improvement, we would expect that the associated neural regions will also show age-related changes. Furthermore, the protracted period of development of the frontal lobes (Giedd et al., 1999;Schmithorst et al., 2002;Giedd, 2004; including connections with Brodmann's areas 17, 18, 31, and 32) may make the associated cognitive functions, the underlying regions, and the connections with these regions particularly dynamic through the course of childhood.
To generate a verb that is plausibly related to the noun, the child must select semantic concepts that are associated with the meaning of the noun. On the output side, semantic retrieval is likely to engage the middle and inferior temporal regions (semantic knowledge) and the hippocampi (information retrieval). The semantic concepts must be coded into phonological form, typically thought of as the second stage of the word retrieval process (Binder et al., 2008). Because semantic associations refine over the course of childhood (McDonald, 1997;McGregor et al., 2002;Beitchman et al., 2008), it is likely that activation in both of these areas will change with age. The phonological form is further coded into subvocal speech (Thompson-Schill et al., 1997). This suggests a second activation by inferior frontal gyrus for subvocal phonological encoding as well as contributions by the insula for speech coordination in subvocal naming. In the covert verb generation task the speech motor network is still engaged but must be inhibited so that words are not spoken overtly (Skipper et al., 2005). We would also expect age-related changes in the neural networks supporting these cognitive components and the connections between them.
There is evidence to suggest that the developing brain adheres to two rudimentary principles of organization: functional integration and functional specialization (Berl et al., 2006). At a neurobiological when compared to the model-based methods that may not identify brain areas with temporal behavior that is not correlated with the experimental design matrix. However, ICA generates a considerable number of components that may not necessarily be part of the studied network (Calhoun et al., 2001b). To address this issue, we have incorporated several additional steps in our ICA method that make the results of our study more targeted and objective . In particular, we have adopted a theory-driven focus based on the Wernicke-Geschwind model of the language network with the aim of investigating developmental shifts in the verb generation circuitry in children from 5 to 18 years of age (Geschwind, 1965a;Anderson et al., 1999). Inclusion of such a focus yields a biologically plausible network model for covert verb generation predicted by the methods proposed here, which is more inclusive and specific in comparison to models extracted using general linear modeling (Yuan et al., 2006;Holland et al., 2007).
In the present paper, we explore the age dependency of the connections between the nodes of the language network that sub-serve the verb generation task. The verb generation model discussed here is based on nodes generated from a group ICA of fMRI data obtained from 336 children ranging in age from 5 to 18 years: which has been discussed in detail in one of our previous publications . The current analysis takes the previously described ICA analysis further by investigating interactions within the identified verb generation network using linear structural equation modeling (LSEM). The initial investigation of verb data using group ICA specifically dealt with: (1) ICA decomposition of verb data; (2) age effects in the task-relatedness of each individual IC map (at individual network level) using an a priori criterion (e.g., correlation with the task reference) and Bayesian formalism; (3) brief description of the steps leading to the present model . In the current analysis, the previously described model is further expanded with input functions (processing of presented nouns) and output functions (retrieval and covert verb production) together with hypothesized connections within and between them. In addition, a theory-driven focus has been proposed taking the biological plausibility of the verb model into account; evaluated against the literature with the focus on the language circuitry. Thus, the current investigation provides an elegant methodology capable of providing unique insight into the framework of neurocognitive brain development in children: by combining and extending previously published ICA results  with some of the figures reproduced in this paper for convenience and completeness. Thus, in the present work the emphasis is placed on the neurocognitive brain development and on estimating the age dependency of inter-and intra-network connectivity predicted using LSEM and correlation analysis.
The verb generation task begins with an auditory presentation of a noun: requires the listener to process the noun's phonological form, and attach meaning to that form (Hickok and Poeppel, 2004). Based on our previous description, this process begins with an input via the superior temporal gyrus . Processing in this area extends posteriorly from primary auditory cortex, as the act of accessing word's meaning is presumed to activate a broader network (Levelt et al., 1999;Pulvermuller, 2001). The posterior superior temporal cortex is known to be structurally asymmetric from the Karunanayaka et al. A developmental model for language Frontiers in Systems Neuroscience www.frontiersin.org generation networks extending the finding of our previous study . The previously described group ICA provides a complete recipe of the prerequisite steps involved in ICA decomposition to identify the key elements underlying a biologically plausible neural network that sub-serve a specific neurocognitive task (Schmithorst and Brown, 2004;Karunanayaka et al., 2010Karunanayaka et al., , 2011Kim et al., 2011). A complete age and sex breakdown of the included subjects (native, monolingual, English speakers) is detailed in Table 1. Based on the Edinburgh Handedness Inventory (Oldfield, 1971), 311 subjects were righthanded, 24 left-handed, and 1 ambidextrous. All subjects were prescreened for any conditions which would prevent an MRI scan from being acquired

FunctIonal ImagIng
All images were acquired using a Bruker 3T Medspec (Bruker Medizintechnik, Karlsruhe, Germany) imaging system. An MRIcompatible audiovisual system was used for presentation of the stimuli. Details of the techniques used to obtain fMRI data from younger children are discussed elsewhere (Byars et al., 2002). EPI-fMRI scan parameters were: TR/TE = 3000/38 ms; 125 kHz; FOV = 25.6 cm × 25.6 cm; matrix = 64 × 64; slice thickness = 5 mm. Twenty-four slices were acquired, covering the entire cerebrum. One hundred ten whole-brain volumes were acquired (with the first 10 being dummy scans) in 5 min 30 s. Techniques detailed elsewhere (Byars et al., 2002) were used to acclimatize the subjects to the MRI procedure and make them comfortable inside the scanner. A whole-brain T1 weighted MP-RAGE scan was also acquired for anatomical co-registration.

Verb generatIon task
The fMRI paradigm of silent verb generation (Holland et al., 2001 is a 30-s on-off block design. All stimuli were presented using MacStim (White Ant Software, Melbourne, VIC, Australia) at a rate of one noun every 5 s, for six stimuli during each 30 s epoch. During the active epochs, the subjects silently generated level, this distinction can be identified with progressive specialization or focalization reflecting the consolidation of synaptic reinforcement of a network (Lenneberg, 1967;Muller et al., 1998;Berl et al., 2006). In this paper, we present a unified framework and examine the developmental trajectories in the language circuitry based on fMRI data using complementary modeling approaches. As previously , we employ ICA, a data-driven method, to identify spatially coherent activation patterns. In the current investigation, we extend these analyses by applying correlation analysis and LSEM to model connectivity between these spatial distributions. Several, other approaches have previously been proposed to investigate network interactions following ICA analyses. Stevens et al.'s (2007) used dynamic causal modeling (DCM) to search for the presence of a meaningful causal structure among selected IC time courses in an event related fMRI study of visual Go/No-Go task. Another study examined the functional network connectivity (FNC) between schizophrenia patients and healthy controls based on the temporal dependency among ICA components (Jafri et al., 2008). Demirci et al. (2009) extended this analysis one step further by incorporating Granger causality test (GCT) to investigate causal relationships between brain activation networks; we also have recently implemented Granger causality analysis to investigate the connections within the epileptic network . Several, other investigations have further highlighted the usefulness of combining ICA with Granger causality on simulated, single subject and group data (Londei et al., 2006(Londei et al., , 2007(Londei et al., , 2010. Some of the above mentioned methods are relatively sophisticated and more suitable for investigating specific group differences between healthy and patient populations. However, the emphasis of the current analysis is on investigating the overall developmental trends associated with the language circuitry and presenting the findings in the framework of a theoretical model for neurocognitive brain development. Thus, given the large sample size, the simplicity of the proposed partially data-driven approach can be considered more suitable for understanding the global network structure (including connectivity) associated with complex verbal language tasks in general, and verbal fluency tasks in particular.

materIals and methods subjects
One hundred sixty-five boys and 171 girls took part in the study following Cincinnati Children's Hospital Institutional Review Board approval. Informed consent was obtained from parent or guardian, an assent was also obtained from subjects 8 years and older. Exclusion criteria were: previous neurological illness; learning disability; head trauma with loss of consciousness; current or past use of psychostimulant medication; pregnancy; birth at 37 weeks gestational age or earlier; or abnormal findings at a routine neurological examination performed by an experienced pediatric neurologist. All subjects were considered healthy based on neurological, psychological, and structural measures . Subjects included in this report were also included in our previous studies focusing on verb generation in children (Holland et al., 2001Karunanayaka et al., 2010) and adults (Szaflarski et al., 2006a). While this report includes fMRI data from the same subjects, it describes an entirely new analysis of connectivity ( causality) associated with covert verb basis for completion of the task through the age of the oldest subjects (18 years). By concatenating the data from the entire cohort and searching for the components (or networks) that are common across the age group, we are able to identify the persistent structural elements (network nodes) underlying the fMRI verb generation task consisting of seven IC networks and shown in Figure 1. Table 2 contains a summary of the respective activation foci for each of the components. Coordinates listed for each IC correspond to the center of mass of each individual spatial element contained in the IC map.
To summarize, the selection of IC maps is based on three criteria: (1) power spectral analysis at the task frequency; (2) phase and (3) relevance of the spatial maps to the theoretical model of verb generation. Thus, by following the above mentioned criteria, the results of this ICA analysis can be replicated by other researchers in the field.
From this point on, we focus on estimating the changes in connectivity between elements of the model identified by ICA using LSEM which is a unique contribution of this work. appropriate verbs such as drink or fill, to aurally presented nouns such as cup. Subjects were asked to tap their fingers in response to a modulated tone presented at 5 s intervals during the control epochs. The control task was specifically designed to control for sublexical auditory processing and also to divert subjects to stop generating verbs into the control epochs. The fMRI task was selected such that children as young as 5 years old would be readily able to perform the task without any difficulty.

group Ica
A complete description of the group ICA methodology for verb generation fMRI data has been discussed in detail elsewhere (McKeown et al., 1998;Calhoun et al., 2001a;Karunanayaka et al., 2010). Basic steps involved in ICA decomposition are briefly mentioned here for the purpose of completeness. ICA is a data-driven analysis technique that does not rely on any prior knowledge of the task performed and is capable of identifying spatially independent components that have similar time courses. The power of group ICA in making statistical inferences from fMRI data has been presented in several investigations (Calhoun et al., 2001a;Schmithorst and Brown, 2004;Karunanayaka et al., 2010).
The ICA decomposition entails several preprocessing steps [normalizing (mean centering) and 40 retained principal components (PCA)] at the single subject level. The data from all subjects are then concatenated into a single dataset before a second PCA reduction resulting in 50 retained components. Finally, 25 runs of the Fast ICA algorithm (Hyvarinen, 1999b) are combined with hierarchical agglomerative clustering (Himberg et al., 2004) to estimate and validate the independent component maps sub-serving covert verb generation. Performing multiple runs (when combined with hierarchical agglomerative clustering) ensures that our analysis resulted in the most reliable components even after taking into account the stochastic nature of the Fast ICA algorithm. Although, ICA can be used to remove motion-related artifacts, individual motion has been fully characterized before performing the ICA decomposition. A detailed analysis of motion (including task-related movement) related to this task is discussed elsewhere (Yuan et al., 2009).
The task-relatedness of each IC map is then investigated using the associated IC time course by examining the spectral power at the task frequency and the phase of the IC time course relative to the task reference function as detailed previously . It should be noted that, by definition, spatial ICA requires independence only in the spatial domain and not in the time domain. Thus, an analysis performed in one domain (e.g., time) can be followed by analysis in another domain (e.g., spatial) without adding any undue bias to subsequent statistical manipulations. Finally, a voxel-wise random effects analysis (one-sample t test) is performed on selected individual IC maps in the spatial domain to determine the cortical regions active in the entire cohort. To further clarify this step, if one were to reverse the domains of the preceding analysis (e.g., spatial followed by time), the end result would be the same because of the above mentioned symmetry. An assumption inherent in this approach is that the structural components of the network for verb generation are in place by the age of the youngest subjects in our cohort (5 years) and continue to get fine tuned to form the structural FIguRe 1 | Seven task-related spatial independent components maps are shown in panels a-g. These ICs are computed using group ICA analysis of 336 children ages 5-18 performing the task of covert verb generation . Slice range: Z = −25 to +50 mm (Talairach coordinates). Three corresponding single subject IC maps are shown at bottom (g, b, d). These individual spatial maps and the associated time courses (Figure 2B) are estimated using a back propagation algorithm following the ICA decomposition at the group level and used in the subsequent LSEM analysis. All images are in radiologic orientation.  . Specifically, representative average time courses were extracted from the functional data set (i.e., real signal intensities) based on these functionally defined ROIs. It is important to remember that ROIs derived from spatial IC maps often include multiple anatomical brain areas, as outlined in Table 2. For all of the IC maps shown in Figure 1, except for IC d, ROIs were defined separately for the left and right hemisphere components of the IC. Based on these ROIs, as mentioned earlier, extracted real signal time courses from the functional data set were then used for the between hemispheres intra-network connectivity computations. A variety of models can be tested in SEM to capture relationships among variables and can provide a quantitative test for a hypothesized theoretical model. SEM takes the entire variance-covariance structure into consideration when evaluating models. Furthermore, SEM is a generalization of regression, path and confirmatory factor models that have been extensively used in psychology, economics and other social sciences. The model estimation in SEM involves minimizing the difference between the observed variancecovariance structure and the one predicted by the implied model. However, when using SEM to model brain activity no distinction is made between the neuronal and the hemodynamic levels (Penny et al., 2004b) which can be considered a drawback of the method.
In the presented model, which is based on Figure 1, we only evaluated the feed-forward connections. As noted above, representative time courses for each of the components (elements) in the LSEM are comprised of individual IC time courses from the previously performed ICA decomposition. The individual LSEM(s) were then solved for optimal path coefficients using the Amos software (Arbuckle, 1989) which utilizes an iterative maximum likelihood method. These optimal path coefficients (connection strengths) correspond to the solution of the structural equations where the difference between the observed and the predicted covariance matrix is a minimum. Finally, we evaluated the goodness of fit between the predicted and the implied covariance matrices using the χ 2 distribution with m (m + 1) − n degrees of freedom (m corresponds to the number of elements and n corresponds to the number of coefficients in the LSEM respectively). The details of LSEM implementation for fMRI data have been discussed elsewhere (McIntosh and Gonzalez-Lima, 1994;Solodkin et al., 2004;Karunanayaka et al., 2007;Dick et al., 2010). The LSEM itself was used (constrained by the proposed verb generation model discussed in the introduction) in a semi-exploratory manner when selecting the final LSEM. Advantages of alternative methods for brain activity modeling (such as DCM) have also been discussed by other authors (Friston et al., 2003;Penny et al., 2004b). Recently, an extended version of SEM called unified structural equation modeling (uSEM; Smith et al., 2010) has been proposed capable of estimating contemporaneous as well as lagged effects simultaneously (Stoeckel et al., 2009). An automatic search procedure has also been proposed to uSEM making it entirely data-driven by increasing its flexibility substantially ). However, DCM is still appears to be the most statistically sophisticated approach that incorporates neuronal hemodynamic relationship into a dynamic model of BOLD activities using Bayesian estimation (Friston et al., 2003;Friston and Stephan, 2007;Sarty, 2007). Thus, given the fact that the relationship between BOLD signal and neuronal activity lInear structural equatIon modelIng Linear Structural Equation Modeling is a statistical method mainly used for hypothesis testing regarding causal influences among measured or latent variables. In addition, SEM is capable of statistically testing a variety of theoretical models that hypothesize how sets of variables define constructs and how these constructs are related to each other. In terms of neuroimaging, SEM relates to effective connectivity that captures causal relationships (directionality) in terms of path coefficients in the model. This approach differs from a typical functional connectivity analysis that can only determine the degree to which two brain regions co-vary (Friston et al., 1997). As mentioned elsewhere, our group ICA decomposition is based on the methods developed by Calhoun et al. (2001b) and is designed to evaluate individual IC maps and corresponding time courses based on group results. In other words, in this method individual IC time courses are estimated using a back propagation method which is followed by the ICA decomposition at the group level. In this paper, we use these individual IC time courses as input to estimate LSEM(s) at the subject level in order to examine the effective connectivity within the network model for verb generation.
A second level, intra-network functional connectivity analysis was also performed using representative real signal intensity average time courses from ROIs defined based on the spatial distribution IC based on known functional neuroanatomy (knowledge-base) ascribed to each Brodmann's area encompassed by the component. As mentioned previously, a theory-driven focus (Geschwind, 1965a;Anderson et al., 1999) complements data-driven methods such as ICA by way of corroborating prior hypotheses about cognitive functions sub-serving the verb generation fMRI task. Depending on the modularity (or function), IC modules are then connected to one another to form the LSEM. In this study, LSEM is directly derived from the covert verb generation model as discussed in the Section "Introduction". For studies of developmental changes within a network, LSEM of an fMRI task can investigate what changes in functional connectivity explain the neural basis of development in language networks. This physiological approach should be guided by the weak constraint that anatomical proximity and connectivity of brain regions are incorporated in the model . Alternatively, a cognitive approach can also be implemented to investigate how functional/effective connectivity changes are related to cognitive development. The emphasis of the current analysis is inline with the latter approach where the effective connectivity changes between IC modules subserving covert verb generation are investigated.
Finally, a second level Pearson correlation analysis was performed on path coefficients in the LSEM to investigate any age effects associated with the proposed cognitive model for covert verb generation.

results
Six out of the seven IC maps shown in Figure 1 were detected in all 25 IC runs while the component shown in Figure 1a was detected in 17 IC runs assuring high reliability  and defines the covert verb generation network for each subject included in the study. The maps in the lower row (individual subject level) of Figure 1 shows three corresponding individual subject level IC maps with corresponding IC time courses: estimated following the ICA decomposition at the group level and used in the subsequent subject level LSEM analysis. Figure 2A shows two of the corresponding average time courses for IC maps shown in Figures 1a,d. Figure 2B shows the individual IC time courses for these networks in two subjects: used as the input to the LSEM computations. The phase progression of the average time courses from leading to lagging the task reference time course (indicated by dark and light gray background) is clearly visualized in Figure 2A.
The developmental trajectories, network behavior (lateralization, task-relatedness, etc.) and the language functions attributed to each IC have been discussed in detail elsewhere . The highly left-lateralized IC map shown in Figure 1d (with lateralization index equal to 1) was identified previously as capturing most of the left-dominance observed in a standard GLM analysis for the covert verb generation task . To perform the intra-network connectivity analysis for this leftlateralized network, four separate ROIs were defined in the left hemisphere as shown and labeled in Figure 3. As explained above, only the real signal time courses from activated regions (refer to Table 2 for further details) inside the colored circles are included in the ROI analysis. The connection between (1) left middle temporal gyrus (LMTG) → (3) left middle inferior frontal gyrus (LMIFG) showed significant age dependent connectivity changes (r = 0.15, is poorly understood (de Marco et al., 2009), LSEM may be a very effective method for making inferences about changes in the causal structure from fMRI time series data.
In addition, several methods have been employed to obtain representative time courses for the components included in a SEM analysis: one popular method being the maximum active voxel representation (Jennings et al., 1998;Goncalves et al., 2001) which we employed previously to investigate developmental trends associated with the narrative story comprehension in children . However, in the current analysis, IC time courses were used to evaluate individual LSEMs to investigate the verb generation task in children. A brief description of the differences between the two methods are included in the section below and also in the Section "Discussion."

bIologIcal constraInts
Several principles have guided the process of constructing a biologically plausible linear structural equation model for verb generation. As the first step (described above), ICA was used as a data-driven descriptor of neural elements involved in performing the fMRI paradigm. The second step involved a Fourier method in the time domain to determine which ICs were most task-related by testing the correlation between the fundamental frequency of each IC time course and the task frequency. The third step involved constructing a biologically plausible LSEM using the knowledge of the sequence of neurocognitive functions involved in the task with IC modules as building blocks . The IC maps require only independence in the spatial domain allowing highly correlated temporal structures to form the theoretical basis for the current SEM analysis. Finally, the phase of the Fourier transform of the associated IC time courses and the known neuroanatomical constraints were also taken into consideration when imposing connections between the model elements.
Some individual ICs out of the seven selected, contain more than one Brodmann's area even though the representative time course for the IC represents all of the voxels included in the spatial map. This is because ICA reveals a set of chronoarchitectonically identified areas (Bartels and Zeki, 2004) or functionally connected regions that may span several Brodmann's areas. If a given cognitive task recruits only one of the observed regions in a given map, then there will be another component separated out by ICA containing only that region. However, if two distinct cognitive functions have very similar time course, they may well be grouped into a single ICA component. This is a limitation of correlational analysis. Still, under certain minimal assumptions, the spatial independence of IC maps can be equated with their modularity, establishing a correspondence between the IC component and a specific cognitive task (Duann et al., 2002;Calhoun et al., 2004). A limitation of this assumption is our inability to determine spatial independence of components with absolute certainty due to the finite number of voxels in fMRI experiments. However, this limitation may have only a minimal effect on the current investigation because of the excellent signal to noise ratio provided by the large number of subjects in the study. Therefore, we argue that it is reasonable to assume that each IC map constitutes a module (a cognitive functional unit) in the proposed LSEM. Depending on the spatial distribution of the IC (activation), a specific language function can be assigned to each Karunanayaka et al. A developmental model for language Frontiers in Systems Neuroscience www.frontiersin.org As described in the Section "Materials and Methods," an LSEM was constructed using the functional IC maps with reference to the literature for prior knowledge (i.e., knowledge-base) about the known neuroanatomy of the brain regions involved in the language circuitry. This LSEM was further refined based on the hypothesized cognitive functions associated with the brain regions encompassed within each spatial IC map, forming the basis for the proposed theoretical cognitive model for covert verb generation as shown in Figure 5. Of note is that the connections between brain regions may not be explicitly included in the proposed model if they are implied by inclusion within a single IC. For example, IC d includes frontal, temporal, and parietal regions. The cartoon in Figure 5 demonstrates this aspect of the model by using an extended ROI spanning these lobes to illustrate the spatial extent of IC d. Table 3 shows the average value of each standardized path coefficient and the age-related changes in path coefficients computed for the LSEM shown in Figure 5. A similar figure (model) was included in a previous study by Karunanayaka et al. (2010); though that diagram did not include the path coefficients computed here as a parameter expressing brain connectivity. As mentioned earlier, the focus of the present analysis is on developmental changes in connectivity within the neural circuitry of language; therefore, we p = 0.007). The Functional connectivity between (1) LMTG → (4) left angular gyrus (LANG) showed no significant age effects. Similarly, the functional connectivity between (4) LANG → (2) left inferior frontal gyrus (LIFG) showed significant age effects (r = 0.143, p = 0.0089) while the functional connectivity between (4) LANG → (3) LMIFG did not. Finally, the functional connectivity between (3) LMIFG → (2) LIFG showed a highly significant age effect (r = 0.18, p = 0.002).
Similarly, we also examined the inter-hemispheric functional connectivity based on individual spatial IC maps. The IC map shown in Figure 1c showed a highly significant age effect (r = −0.3, p = 2.457e − 008) in the connectivity between the hemispheres (Figure 4). Similarly, the IC map shown in Figure 1f also showed significant age effect (r = −0.132, p = 0.015) in inter-hemispheric connectivity. However, the IC shown in Figure 1e (bilateral superior temporal gyri; BA 22) did not exhibit significant age effects in functional connectivity between the left and right hemispheres (Figure 4). Similarly, IC modules a, b, and g also did not exhibit any age-related inter-hemispheric functional connectivity changes. Thus, for these components, we have not included the results of the above mentioned inter-hemispheric functional connectivity analysis.  attributed to the youngest subjects having higher than average IQ . To be more specific, when the children between the ages 5 and 8 years were excluded, the weak correlation between age an IQ did not reach significance and consequently we have not included IQ as a covariate in the analysis.

dIscussIon
Methods for network connectivity analysis based on functional neuroimaging data are developing rapidly as a means of expanding our understanding of neurocognitive function beyond what the neo-phrenology or functional blobology of fMRI have been able to reveal (Friston et al., 2003;Rajapakse et al., 2008;Dick et al., 2010). ICA is an ideal preliminary step for network connectivity analysis because it is able to detect areas that exhibit task-related behavior which might not correlate highly with an a priori model or reference function. In the present analysis, we began with ICA of verb generation data which detected activations in multiple networks with different temporal signatures. Multiple activation time courses detected in the same brain regions (specifically frontal and temporal cortex) provide direct evidence of their participation in multiple cognitive aspects of the verb generation task. ICA provided the basis for construction of a LSEM for the network that sub-serves verb generation task and allowed us to use this standard statistical methodology to explore the age dependency of the relationships among cognitive modules revealed by the ICA analysis . The theoretical framework guiding this research focuses on investigating the developing brain from a network perspective and lays the foundation for deciphering any developmental trends as interactions between underlying networks. Starting with the Wernicke-Geschwind model for the language network, we used a data-driven approach to analyze results from an fMRI experiment in a large sample of children over a wide age range in order to extract key network elements supporting verb generation. This classical model guided our thinking about how to connect modules identified by group ICA results as having a strong correlation with the task behavior. We then examined the network structure to identify developmental trajectories that correlate with age and ability of children to think and reason at increasing levels of maturity . We have shown elsewhere, how ICA can be used to explore developmental changes in brain activation patterns associated with individual neural networks supporting covert verb generation . The current analysis takes this approach one step further by incorporating LSEM to the investigations of the theories of brain development using the regionally weighted or focal network models (Berl et al., 2006). Although, these hypothesized brain developmental models draw support from current neuroimaging literature, our analysis seems to favor the regionally weighted model of normal language development.
Independent component analysis by itself is not capable of revealing the precise cognitive correlates of the identified components . Instead, this data-driven method must be utilized to identify spatial distributions (IC maps) from fMRI data. As with GLM-based analyses, the function of the detected regions must be inferred and should be constrained by prior knowledge of examined changes in the path coefficients estimated by the LSEM as a function of age. The following path coefficients exhibited agerelated changes: The path coefficient between IC e → IC f showed an increase in connectivity with age (r = 0.13, p < 0.017). The path coefficient between IC e → IC d showed a modest (identified with a trend) age-related connectivity decrease (r = −0.111, p < 0.044). However, the path coefficient between IC f → IC d exhibited a highly significant age-related increase in connectivity (r = 0.18, p < 0.00088). Figure 6 graphically displays the corresponding standardized path coefficients that showed statistically significant age-related changes. These values are italicized in Table 3.
For the group of children included in this study, subject age was significantly correlated with the full-scale IQ (Spearman's r = −0.18, p < 0.0008). This small but significant negative correlation is mainly

FIguRe 3 | Regions included (only the active areas) in the intracomponent functional connectivity analysis for IC d are 1 medial temporal gyrus (LMTg), 2 inferior frontal gyrus (LIFg), 3 middle inferior frontal gyrus (LMIFg), 4 angular gyrus (LANg).
Each brain region will be represented by the average activation within that ROI across time. Slice range: Z = −25 to +50 mm (Talairach coordinates). All images are in radiologic orientation (left in the picture is right in the brain). Figures 1C,e. IC c exhibits a highly significant functional connectivity between the left and the right inferior frontal gyrus.

Karunanayaka et al. A developmental model for language
Frontiers in Systems Neuroscience www.frontiersin.org . Given the limitations [(Wright's rules; Write, 1934) and the number of nodes in the model] in evaluating LSEM(s), a careful consideration must be given before selecting either approach. In general, any theoretical model for language related cognitive functions will be a compromise between the complexity of the neural system sub-serving language comprehension and the interpretability of the resulting models. Complex models can account for intricate dependencies (both anatomical and functional) but the interpretability of the resulting models would be severely compromised (McIntosh and Gonzalez-Lima, 1994;Dick et al., 2010). As suggested by Dick et al. (2010), one approach would be to use the hypotheses being tested as guiding the constraining aspects of the model development. An alternative, more appealing approach would be to model brain functions in terms of interactions between underlying sub-networks, inline with the method we have proposed in this paper. To circumvent inherent drawbacks of the second approach, in addition to the theory-driven focus, we incorporated a secondary correlation analysis specifically to investigate the within network behavior sub-serving covert verb generation in children (Friston et al., 1997).
The functional connectivity results of IC d revealed unique features related to semantic processing circuitry in children. Several studies have implicated activation in the middle temporal gyrus the functional neuroanatomy. However, once the spatial distributions are known, depending on the complexity either a physiological or a cognitive approach can be employed for the connectivity analysis  ; Figure 4). The brain cartoon shows the approximate locations of each IC map from Figure 1. Transparent ellipses indicate regions located medially within the brain and not visible from the lateral surface whereas opaque ellipses correspond to regions that are mainly located on the lateral surface of the brain. IC d is represented in both frontal and temporal-parietal regions as reflected in the distributed nature of this left-lateralized network. The network is divided into word processing (shown in blue) and word generation modules (shown in green). The SEM block diagram at bottom shows how these brain networks are graphically connected forming the basis for the cognitive model for the covert verb generation task. Only the Feed Forward Connections are evaluated.   lateralization and localization over the course of language development (Ahmad et al., 2003;Gaillard et al., 2003). Our previous findings of increasing left lateralization of IFG activation with age for the verb generation task in children are consistent with the functional connectivity findings showing decreasing left-right connectivity with age suggesting that the left hemisphere is able to act more autonomously in support of word generation as the brain matures . This interpretation is also consistent with the regionally weighted model of normal language development. Further, this finding alone can explain the differences between young and old subjects in language recovery after left-hemispheric injury with the ability of the language functions to shift to the right hemisphere in the early (prenatal and early postnatal injury) but dependence on the left-hemispheric regions for aphasia recovery in late life stroke (Tillema et al., 2008;Saur et al., 2010). The inter-hemispheric functional connectivity between the posterior aspects of superior temporal gyrus (IC e) showed no age effects. The time courses for IC e and IC c described above have shown the highest increase in task-relatedness (developmental trend) as detailed in a previous study involving the same subject population . However, the age dependence of these networks differs in terms of the inter-hemispheric connectivity as seen in Figure 4, with no significant age trend found in the posterior network encompassed by IC e (BA22).
Although the relationship between structural maturation and functional activation is rather complex, the present functional connectivity data provides additional evidence in support of language lateralization being dominated by the inferior frontal brain regions. While one recent study did not observe any asymmetries in language lateralization in newborns (Kotilahti et al., 2010), this study also found a more uniform involvement of the left hemisphere in speech processing indicating that left-hemispheric specialization for language processing may already be present at birth. Dehaene-Lambertz et al. (2002) also found that left lateralization of language function was present in posterior brain regions in infants as young as 3 months of age. These findings are inline with previously reported left lateralization of language functions noted in 6-to in the acquisition of semantic representations (Blumenfeld et al., 2006;Booth et al., 2007). Similarly, research in adults suggests that more activation in the inferior frontal cortex is associated with more effortful retrieval or greater selection demands (Seger et al., 2000;Gurd et al., 2002;Whatmough et al., 2002;Booth et al., 2007). Age effects seen in the functional connectivity between these two regions suggest that the selection demands imposed on the inferior frontal gyrus increase with age. This may be due to the fact that the present verb generation task does not impose restriction on the number of verbs a subject can generate for a given noun. Evidence suggests that this design is successful in minimizing the amount of variance attributable to performance (Gaillard et al., 2003).
The functional connectivity between LMTG → LANG showed no significant age effects. The inferior parietal cortex has been implicated in feature integration and semantic categorization to form a coherent concept so that semantic relationships between words can be determined (Grossman et al., 2003;Karunanayaka et al., 2010). The demand for such processes may be at a minimum for this task (ceiling effect) since we developed this fMRI task in such a manner that even the youngest children in our study can perform this task easily. Nevertheless, Booth et al. (2007) have suggested that the inferior parietal lobule may have distinct areas for processing semantic versus phonological information. This may explain observed age effects in functional connectivity: between LANG → LMIFG with no age effects and between LMIFG → LIFG with highly significant age effects.
The significant decrease in functional connectivity with age between right and left hemisphere elements of IC c implies a substantial change in the degree to which the left and right brain regions (inferior frontal gyrus) co-vary. Note that structural and functional asymmetries have also been found in the prenatal and early postnatal brain (Wada et al., 1975;Chi et al., 1977;Dehaene-Lambertz et al., 2002) suggesting a bias for left hemisphere language lateralization very early in life. The anatomical data suggest that early brain development may lead to an underlying architecture that preferentially supports language within the left hemisphere: a normal variant of the focal network model (Berl et al., 2006). This neuroanatomical bias is hypothesized to be related to functional  verb generation in the developing brain. As mentioned previously, compared to DCM, the current investigation only models contemporaneous connections without taking into account the neuronal hemodynamic relationships explicitly (Penny et al., 2004b). The emphasis is, therefore, on the overall network behavior confirming or facilitating the generation of new hypothesis. The current investigation focused on the overall connectivity pattern shedding more insight into several networks that need further investigations using more sophisticated methods like uSEM, DCM, or Granger causality (Stevens et al., 2007;Jafri et al., 2008;Demirci et al., 2009;Londei et al., 2010;Smith et al., 2010). SEM is useful in this regard in that it provides a quantitative measure of overall model fit which allows the optimum set of path coefficients to be identified objectively. These coefficients can then be examined as a function of age to determine how connection strengths change with brain development. Finally, LSEM was also used as an exploratory tool in the proposed theoretical model in a highly restrictive manner. By introducing a theory-driven focus we partially avoided evaluating models of different structures. However, model selection (or identifying the true network structure) is a challenging statistical problem that has received increased attention in the neuroimaging community in recent times (Zheng and Rajapakse, 2006;Rajapakse and Zhou, 2007). We have already developed a Spectral Bayesian Network method (based on Model Averaging) to identify the most plausible models based on fMRI data, which is inline with our long-term objective of developing statistical methods capable of confirming (or rejecting) existing theoretical models for cognitive development.
lImItatIons Study limitations inherent in covert verb generation task have been discussed in detail elsewhere (Szaflarski et al., 2006a,b;Karunanayaka et al., 2010). Therefore, we will only review additional limitations pertaining to the analyses employed in this paper.
In this study, we have only focused our attention on task-related networks even though considerable amount of intrinsic fluctuations are typically inherent in fMRI time courses. ICA, in general, tends to over specify the problem imposing severe limitations on our computational ability for connectivity analysis. Implementing objective methods to select non-task-related components to be included in the connectivity analysis is non-trivial. On the other hand, including such components is very subjective making interpretations difficult. The SEM should also be limited to a reasonable number of nodes (maximum of 10∼15) as any data set can be fitted to models with increasing complexity. Thus, we have adopted a theory-driven focus coupled with proper selection processes to guide the analysis and interpretations circumventing above mentioned drawbacks. Therefore, we had no option but to limit the analysis to task-related components. However, if one can overcome the computational (and methodological) limitations, DCM might be more suitable to investigate intrinsic connectivity that is affected by the context of the task in ways which do not show up as a strictly task-related modulation of the time course.
As previously mentioned, ICA is a data-driven technique and, therefore, its use obviates conventional statistical approaches to hypothesis testing. Consequently, one extension of this data analysis method would be to incorporate constraints at the ICA 12-month-old children (Minagawa-Kawai et al., 2007) and later studies of language lateralization in older children, adolescents and young adults (Holland et al., 2001;Szaflarski et al., 2006a).
The proposed LSEM for verb generation is hypothesized to support both word processing and word generation. However, only the networks included in the word processing module exhibited age dependent effective connectivity changes. Each of these networks represents a unique spatial distribution with corresponding time course that sub-serves specific functions of the network (e.g., working memory, visual imagery, or acoustic word recognition). As mentioned earlier, although the spatial distributions of IC maps are independent, the corresponding time courses are allowed to have highly correlated temporal structures.
According to the focal network theory, the underlying neural network structure for language processing is generally well established by the age of 5 (Ahmad et al., 2003) with first evidence of network structure seen already in newborns (Kotilahti et al., 2010). Therefore, it is reasonable to assume that interactions between functions such as coordination of speech articulation, subvocal word production, and visual imagery at network level are well established for this group of children. However, based on our results, there is ample evidence to suggest that the within network (intra-network) behavior is undergoing a continuous process of dynamic change. As discussed in detail elsewhere , the areas of a distributed network can change the degree of engagement making it a more efficient component of the normally developing network. This forms the basis for the regionally weighted model and the differences in weights may account for the observed normal variations in cognitive skill level, use of different cognitive strategies and changes in the biological substrate for a function (Berl et al., 2006). This picture is consistent with the intra-component functional connectivity results observed for IC c, d, and e. .
As mentioned above, module IC d is the most left-lateralized part of the network for this task and is presumed to be associated with semantic representations of the nouns that are being heard . All connections to this module are age dependent. This module may also sub-serve working memory required by the verb generation task. Several studies have reported age-dependent BOLD signal and connectivity changes mainly in the frontal areas of the brain (Gaillard et al., 2000;Schlaggar et al., 2002;Schmithorst et al., 2002;Schapiro et al., 2004). We suggest that these later aspects of development are captured by the observed connectivity changes within the word processing module in our proposed model. Finally, in terms of the regionally weighted model, these changes can be interpreted as increasing the participation of this left-lateralized network supporting phonological and semantic expressive functions as part of covert verb generation.
The biological relevance of the model derives from two sources. First, the highly task-related elements of the model are selected based on the data-driven ICA results. Secondly, the biological plausibility originated with the close correspondence to the Wernicke-Geschwind model and has been evaluated against the literature for the neural circuitry of language; especially for the semantic processing network . These biological underpinnings for our model give us confidence that the proposed model is indeed relevant to the cognitive and biological processes taking place during dependent changes in network connectivity. Connections between networks such as word recognition, phonological working memory and semantic processing exhibited age dependent connectivity changes to varying degrees. We also highlighted the importance of combining network analysis at both the local and global levels.
The intra-network connectivity analysis of the highly left-lateralized network (IC d) provided convincing evidence for ongoing within network connectivity changes. These changes highlight the importance of the degree of engagement of each brain region within a network and can be considered a part of normal development. The results show the capacity of ICA to detect additional task-related spatial modules and their developmental trajectories. Age dependence found in both the path coefficients and in the task-relatedness of this network support the idea that this network will be used increasingly for covert verb generation as children mature. Combined with the results of our previous investigations, the current analysis further endorses the regionally weighted model for neurocognitive development that highlights the importance of interactions between brain regions within a network.
The main emphasis of the current analysis was on investigating the age dependent connectivity trends in language networks in the developing brain using functional neuroimaging data and an inferential statistical method capable of testing hypothesis based upon existing theoretical models for cognitive functions. Modeling cognitive functions will always be a compromise between the complexity of the underlying neural systems and the interpretability of proposed models. However, by carefully selecting appropriate analysis methods in an incremental manner of complexity, theories of cognitive development can be investigated in detail. This new approach to network modeling reveals a more complex developmental trajectory than proposed by models currently in use, e.g., the Wernicke-Lichtheim-Geschwind theory (Geschwind, 1965a,b) and has the potential to lead to a deeper understanding of developmental changes in cognitive function.

acknowledgments
The study was supported by a grant from the U.S. National Institute of Child Health and Development, R01-HD38578. The authors acknowledge the assistance of Dr. Anna Byars, Ph.D., in the administration of the Wechsler IQ tests, and Drs. Richard Strawsburg, M.D. and Mark Schapiro, M.D., for performing the neurological examinations. Dr. Szaflarski is currently supported by NINDS K23 NS052468. Thanks to Meghan O'Connell for applying her medical illustration skills to Figure 5. decomposition step to guide the analysis and increase the predictive power (Lu and Rajapakse, 2005). Finally, the verb generation task utilized here was not specifically designed to acquire in-scanner performance data. Consequently, performance effects on the connectivity coefficients cannot be completely discounted even though the fMRI task design enabled the youngest children in the study to complete the task without any difficulty. Since, performance can be related to IQ, including IQ as a covariate can produce overcorrected, anomalous, and counterintuitive findings about neurocognitive functions (Dennis et al., 2009). It has also been shown that IQ should only be used as a covariate in those rare circumstances where selection bias has produced problems of non-representativeness in the sample (Dennis et al., 2009). Clearly, such a condition was not present here although we observed a small negative correlation between age and IQ. This was mainly due to our youngest subjects having higher than average IQ scores. Furthermore, in one of our previous connectivity studies (narrative story comprehension) with the same population, the effects of the age × IQ interaction term were investigated using a multivariate regression model and were found not to confound the age-related tendencies associated with SEM path coefficients . This performance-related limitation can be addressed in the future by collecting intra-scanner performance data using either sparse fMRI data collection  or block-design task with forced responses (Szaflarski et al., 2002). Such a design will also allow real-time performance on the task to be monitored and potentially included as a covariate in the analysis of age dependence in connectivity. Recently we have shown that brain activation during covert verb generation correlates with the number of verbs generated during an overt phase of verb generation during the same task . While both overt and covert verb generation produced similar patterns of activation, the correlation with performance suggests that performance could also be related to connectivity in the language networks sub-serving the tasks. This question could be specifically addressed with a modified overt verb generation task in which the number of responses is explicitly controlled as a design parameter. conclusIon A theoretical model for covert verb generation was investigated using fMRI data from a large cohort of children and adolescents between the ages 5-18 years undergoing fMRI study with such a task. Previously identified, spatially independent and task-related networks (IC maps) were combined with SEM to investigate age