Ahmed Yureidini

Erwan Kerrien

Jérémie Dequidt

Christian Duriez

S Cotin

Local implicit modeling of blood vessels for interactive simulation

published or not. The documents may come

Introduction

In the context of interventional radiology, interest is growing in the potential benefits of computer-based simulators. Overall, the expected benefit is a reduction of time: for training, operation and hospitalization; but also to promote new surgical devices and techniques. Such expectations are heightened in the intravascular treatment of brain aneurysms for which highly skilled practitioners are required. Recent works [START_REF] Dequidt | Towards interactive planning of coil embolization in brain aneurysms[END_REF] attest the feasibility of real-time simulators for such procedures. However, a key requirement is the availability of the blood vessel surface. We propose in this paper an original blood vessel modeling algorithm, specifically suited to real-time simulation of the interventional gesture.

Vessel lumen segmentation can rely on a vast diversity of models [START_REF] Lesage | A review of 3D vessel lumen segmentation techniques: Models, features and extraction schemes[END_REF]. Here, a catheter, or guide, is pushed through the vessels, interacting with their surface. The vessel surface model should enable a smooth motion of the tools and a precise and efficient collision detection. Implicit surface representations, where the surface is defined as the zero-level set of a known function f , are arguably well suited. First, C1-continuous models (with C0-continuous normal) allows for much smoother sliding contacts. Unwanted friction may occur with polyhedral surfaces or level-sets defined on a discrete grid. Second, implicit surfaces offer an improved collision management over parametric surfaces. Indeed, the implicit function value at a point tells whether this point is inside or outside the surface, detecting a collision in the latter case. Furthermore, the implicit function gradient gives a natural direction for the contact force used to handle this collision.

Implicit vessel modeling has been studied, in particular in the context of visualization [START_REF] Preim | 3d visualization of vasculature: An overview[END_REF]. A major problem is to handle unwanted blending between two structures that are geometrically close, but topologically far from each other. Locally adapting the blending [START_REF] Bernhardt | Implicit blending revisited[END_REF] has proven efficient, but is not a valid answer in the context of simulation because the simulated tool motion is discrete. The inclusion test might therefore remain successful even though a crossing over a small gap occurred during the time step. Since the seminal work of Masutani [START_REF] Masutani | Region-growing based feature extraction algorithm for tree-like objects[END_REF] the need for local models, based on the centerline of the vascular tree, has been recognized. Placing a graphics primitive at each point of the centerline [START_REF] Tyrrell | Robust 3-D modeling of vasculature imagery using superellipsoids[END_REF], or using convolution [START_REF] Preim | 3d visualization of vasculature: An overview[END_REF] lacks precision where the vessel section is not circular (or elliptic), and does not correctly handle pathologies such as aneurysms. Schumann [START_REF] Schumann | Implicit vessel surface reconstruction for visualization and CFD simulation[END_REF] addressed both problems and used Multi-level Partition of Unity (MPU) implicits, to get a locally defined model with large and general modeling capabilities. However, the blending issue is not addressed and the model is not related to the topology. Furthermore, MPU implicit gradient gives an appropriate contact direction close to the surface but could mislead contact forces elsewhere.

The proposed approach reconstructs the vessel surface as a tree of local implicit models: one local implicit model is placed at each point on the vessel centerline. Complex branching shapes such as the brain vasculature can thus be modeled in a way that is very efficient for interactive simulation. Implicit surfaces generated by skeletons were chosen for their locality and their genericity. Model fitting is defined as an energy minimization problem, alternating with a model refinement as in [START_REF] Tsingos | Implicit surfaces for semi-automatic medical organ reconstruction[END_REF], but improving over this latter work. Our complete modeling algorithm is described in Section 2. Section 3 reports our results on a set of 5 3D rotational angiography (3DRA) patient data. Finally, we discuss our results and give a conclusion of our work in Section 4.

Modeling algorithm

Our algorithm requires that the vessel centerline is extracted as a tree, not necessarily dense, and that a local vessel radius estimate is available at each point on this centerline. Many valuable methods can provide this initial structure [START_REF] Lesage | A review of 3D vessel lumen segmentation techniques: Models, features and extraction schemes[END_REF]. Moreover, a set of points, located on the local vessel surface, should be associated with each point on the centerline. In practice, the results presented in Section 3 rely on [START_REF] Yureidini | Robust RANSAC-based blood vessel segmentation[END_REF]. However, should the vessel surface be available as a mesh, selecting the vertices in the vicinity of each point on the centerline should also provide such a point set. In this case, our method could be seen as a mesh implicitization. Our overall idea is to fit each local point set with an implicit surface generated by a point-set skeleton (Section 2.1). The local models are organized under the same tree structure as the point centerline (Section 2.2). Solving contact constraints at each point of the interventional tool, is performed using the appropriate local implicit model as the vessel surface (Section 2.3).

Local Implicit Modeling

Implicit formulation An implicit iso-surface generated by a point-set skeleton is expressed as the zero-level set of a function f , a sum of implicit spheres:

f (X; p) = T - N b j=1 α j φ |X-Cj | ρj .
T is the isosurface threshold, {α j } are positive weights, and {C j } is the point set skeleton. Each implicit sphere #j is defined by a symmetric spherical function, centered on C j , of width ρ j . The local field function, or kernel, is a function φ : R → R + , rapidly decreasing to 0 at infinity. For example, all results of Section 3 were produced using the 'Cauchy' kernel [START_REF] Sherstyuk | Kernel functions in convolution surfaces: A comparative analysis[END_REF]:

φ(x) = (1 + x 2 /5) -2
(dividing factor 5 normalizes the kernel such that φ (1) = 0). Such objects were called differently depending on the kernel used [START_REF] Sherstyuk | Kernel functions in convolution surfaces: A comparative analysis[END_REF]. Our method is not kernel-dependent, and was successfully used with the computationally less efficient Gaussian kernel. Muraki [START_REF] Muraki | Volumetric shape description of range data using Blobby Model[END_REF] was the first to use this type of model in the context of object reconstruction. Following this seminal work, we shall use the terms blob for an implicit sphere, and blobby models as a generic name for the implicit models.

In our particular simulation context, in order to help predict collisions, and have the function give a valid contact force direction, the algebraic value f (X; p) at point X should relate monotonously to the geometric distance of X to the surface. We set α j = ρ j , which obviously establishes the sought relation in the case of a single blob. Meanwhile, redundancy in the parameters of f is dismissed.

Energy Formulation Fitting a surface to N points {P i } 1≤i≤Np can be written as an energy minimization problem [START_REF] Tsingos | Implicit surfaces for semi-automatic medical organ reconstruction[END_REF][START_REF] Muraki | Volumetric shape description of range data using Blobby Model[END_REF]. We propose to combine 3 energy terms: E = E d + αE c + βE a where (α, β) ∈ R + 2 , and:

-E d = 1/N p i f (P i ; p) 2
translates the algebraic relation between data points and the zero-level set. It gives a raw expression of the approximation problem.

-

E c = 1/ (N b (N b -1)) j =k s √ ρj ρ k |C k -Cj | 12 -2 s √ ρj ρ k |C k -Cj | 6
is Lennard-Jones energy. Each term is minimal (with value -1) for |C j -C k | = s √ ρ j ρ k , being repulsive for blobs closer than this distance, and attractive for blobs further away. It imposes some cohesion between neighboring blobs to avoid leakage where data points are missing, while preventing blobs from accumulating within the model.

-E a = 1/N p i κ(P i ) 2
κ(P ) is the mean curvature. It can be computed in a closed form at any point in space from the implicit formulation [START_REF] Goldman | Curvature formulas for implicit curves and surfaces[END_REF]

κ(P ) = ∇f t H f ∇f -|∇f | 2 trace(H f ) 2|∇f | 3
where ∇f is the implicit function gradient and H f its Hessian matrix, both computed at point P . This energy smoothes the surface according to the minimal area criterion. In particular, the wavy effect that could stem from modeling a tubular shape with implicit spheres, is reduced.

Behind the rather classical form given above for the energy terms, it is important to notice that the whole energy is known under a closed-form expression. As a consequence, closed-form expressions were derived for its gradients with respect to the blobby model paramaters {ρ j } and {C j }.

Selection-Subdivision

The blob subdivision procedure proposed in the seminal work [START_REF] Muraki | Volumetric shape description of range data using Blobby Model[END_REF] was exhaustive and time consuming. A blob selection mechanism was added in [START_REF] Tsingos | Implicit surfaces for semi-automatic medical organ reconstruction[END_REF], measuring the contribution of each blob to E d in a user-defined window, and choosing the main contributor. User input is not an option in our context where thousands of blobby models are handled (see Section 3. Moreover, we experimentally noted that this technique was prone to favor small blobs, thus focusing on details, before dealing with areas roughly approximated by one large blob. This behavior is caused by this selection mechanism using the algebraic distance to the implicit surface. Our criterion relies upon the geometric distance approximation proposed by [START_REF] Taubin | Estimation of planar curves, surfaces, and nonplanar space curves defined by implicit equations with applications to edge and range image segmentation[END_REF]. Point P i * farthest to the surface is such that:

i * = arg max 1≤i≤Np |f (P i ; p)| |∇f (P i ; p)| (1)
The blob #j * whose isosurface is the closest to P i * is selected (according to Taubin's distance). Note that this criterion is valid in large areas because we set α j = ρ j in the definition of f . The subdivision step then replaces this blob with two new ones. Their width ρ j * is chosen such that two blobs, centered on C j * , of width ρ j * would have the same isosurface as one blob centered on C j * , with width ρ j * (the formula depends on the kernel). The first new blob is centered on C j * , while the second is translated by ρ j * /10 towards P i * .

Optimization Such a gradual subdivision procedure may lead to a dramatic increase in the number of blobs, and hence the size of the optimization problem. The locality of the kernel φ allows us to focus the optimization onto the newly created pair of blobs. More exactly, only the new blob that is slightly misplaced is optimized, the other blobs remaining constant. The energy is minimized using Polak-Ribiere conjugate gradient (PR) algorithm, taking advantage of the closedform expressions of both the energy and its gradients. A single minimization loop consists in one PR minimization over the center (3 variables), followed by one on the width (1 variable). In practice, a maximum of 5 loops proved sufficient.

Global Implicit Modeling

The algorithm is initialized by placing two blobs at each node on the centerline that form an elongated shape along the vessel direction, with the same diameter as the vessel width. The local data point sets are concatenated so as to overlap with their neighbors, and generate smooth transition between adjacent blobby models. The following three steps are then applied on each blobby model:

1. A first energy minimization is performed over the full blobby model, with a single minimization loop (all centers then all widths of its blobs). 2. The subdivision process is applied, as described above. The process is stopped when a maximum number N s of subdivisions is reached, or the distance between P i * and blob #j * drops below a threshold t g . 3. The model is fine tuned by a single energy minimization loop over the full blobby model, as in step 1. Redundant blobs are either ejected far away from the surface during the optimization over the centers, or their width is reduced to almost zero during the minimization over the widths. A simple clean-up procedure was applied after each minimization loop: blobs further than 20 mm from the node (twice the diameter of the largest artery), and blobs whose width was below 0.02 mm (10% of the voxel size) were removed from the blobby model. In combination with the subdivision process, this clean-up enables the algorithm to simulate large blob displacements, which is hardly possible with the local minimization, thereby adding robustness to poor initialization. Figure 1 illustrates the essential steps of this fitting algorithm and its capacity to model complex shapes, such as aneurysms, even from rough initialization.

Note that in this algorithm, each local blobby model is treated independently from the others, enabling a parallel computation.

Using the local models for simulation

Most of the interventional tools are slender and their motion can be defined by that of longitudinal nodes. During a simulation step, contact with the vessel surface is detected and solved for each point on the tool. Each tool point is thus linked to its closest point on the centerline and the associated local blobby model is used as the surface constraint. In order to take into account the topology, only the neighbors of the current centerline point are considered as candidates to update the surface constraint during the motion of the tool (see Fig. 2).

Our model is adapted to simulation, but not to visualization. We translated the above model selection procedure to give a visual impression of the result in the next section: each local isosurface was first extracted with marching cubes, and was then cut by the median planes separating the current centerline point from each of its neighbors. No blending was performed between adjacent blobby models. The final surface presents as a stack of individual surfaces.

Results

Experimental setup

A set of 5 patient data was used for validation. Each patient data set consisted of a 3DRA acquired on a vascular C-arm (Innova 4100, GE Healthcare) during the intra-arterial injection of the internal carotid artery. Each 3DRA volume is a 512 3 isotropic voxel cube, between 0.18-0.22 mm voxel size. For each patient, the vessels were tracked and points were extracted at their surface using [START_REF] Yureidini | Robust RANSAC-based blood vessel segmentation[END_REF].

The implicit modeling procedure was applied to all 5 datasets. The algorithm parameters were tuned on the carotid artery of one patient, and thereafter used for all patients. The 'Cauchy' kernel was chosen for its computational efficiency. The energy weights were tuned so the energy values are of the same order of magnitude, leading to α = 10 -5 and β = 10 -3 . For Lennard-Jones energy, s was set to 2 (natural distance between two blobs is twice the geometric mean of their width). The isolevel value T was set to 0.1. Indeed, since we set the weight of each blob equal to its width, blobs whose width is below T do not generate any isosurface (maximum function value is ρ j < T ). We therefore chose T to be approximately half the voxel size.

Geometric precision assessment

The geometric error of fit of a local blobby model with respect to a point set was estimated using Taubin's approximate distance. When measured at each data point, it provides a set of statistical samples for this error: the error of fit for a blobby model, d bm was defined as the 90th percentile on this set, which is a more strict but also more robust measure than the mean. Table 1 reports the distribution, in percents, of d bm measured on all 42840 blobby models. Four classes were considered: below 0.5, 1, 2 and over twice the voxel size (which varies between 0.18 and 0.22mm). Two parameters in particular have an impact on the geometric precision: the distance threshold t g , expressing a targetted accuracy of fit, and the maximum number N s of subdivision, expressing the maximum complexity allowed for the models. Measurements were performed for t g = 0.5, 0.3 and 0.2 mm, N s = 100, 50, 30 giving 9 algorithm configurations.

At most 1% of blobby models had an error of fit above 1 voxel. t g appeared as the only parameter to impact the result and t g = 0.3mm was the only case where N s had a noticeable influence. In our opinion, N s only has an influence when t g is adapted to the details to be modeled in the data: t g = 0.5mm was too rough and t g = 0.2mm led the blobby models to fit noisy data points. t g = 0.3mm and N s = 100 provided the best algorithm configuration.

Model efficiency and computation time

The computation time to model one patient was between 15 and 30 minutes, depending on the number of nodes in the vascular tree (between 5461 and 13956 nodes). A sample result is shown on Figure 3, demonstrating a very good vessel resolution and smooth transitions between neighboring blobby models, except at some very isolated locations. Discontinuities between models do not depend on the kernel, but are mostly sensitive to inhomogeneities in the point set density. Despite discontinuities are below voxel size, we believe they prohibit using our model for blood flow simulation using Computational Fluid Dynamics (CFD). Therefore, direct reconstruction from image data is currently being investigated. A rule of thumb to measure the efficiency of our model, is to count the number of primitives used. For each patient we compared the total number of blobs (with t g = 0.3 and N s = 100) to the number of triangles in a vessel isosurface mesh of a similar visual quality (see e.g. Fig. 3, left). The ratio #triangles/#blobs ranged from 7 to 13.8, with an average of 10.5. On a more local basis, the average number of blobs per local blobby model was 4.2 (min=3.7, max=4.6).

We have presented a patient-specific blood vessel surface reconstruction method that is particularly suited to the simulation of interventional radiology procedures. The final model presents as a tree of implicit surfaces generated by local skeletons. Each surface has a limited area of influence which overlaps with its topological neighbors to provide a smooth transition during the tool guidance. The model fitting is guided by the minimization of an energy, with a closed-form expression, alternating with an efficient blob selection and subdivision scheme. The reconstructed models provided a sub-voxel approximation of the vessel surface on 5 patient data. This model was integrated in our simulation software and preliminary tests proved the adequacy and efficiency of our model in that context. As a perspective, we aim at reconstructing directly from image data.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Implicit modeling of an aneurysm. The points {Pi} are in red. (From left to right) Initialization with a single blob ; after the first minimization ; after 25 subdivisions ; final result (100 subdivisions)

Fig. 3 .

 3 Fig. 3. Visible assessment (from left to right): isosurface from the raw data set ; isosurface from the blobby models (13956 models, 51272 blobs, tg = 0.3, Ns = 100): much more vessels are visible; close-up showing smooth transitions between models ; close-up on a small artery branching onto the carotid: the connection is difficult to model

Table 1 .

 1 Distribution in % of the blobby models according to the error of fit (d bm , in voxels) in 4 classes (left column). 9 algorithm configurations were investigated depending on parameters tg (targetted accuracy of fit) and Ns (maximum model complexity).

	tg		0.5mm			0.3mm			0.2mm
	Ns	30	50	100	30	50	100	30	50	100
	d bm < 0.5	47.44 47.18 47.29 60.28 65.07 70.52 89.86 90.07 90.14
	0.5 < d bm < 1 51.64 51.82 51.83 39.17 34.53 29.11	9.75	9.60	9.57
	1 < d bm < 1.5	0.85	0.96	0.83	0.48	0.34	0.32	0.32	0.26	0.23
	2 < d bm	0.06	0.04	0.04	0.06	0.06	0.05	0.08	0.07	0.06

Acknowledgments

The authors are grateful to Prof. R. Anxionnat (CHU Nancy) for the data.


Bruno Chanteau
email: bruno.chanteau@univ-paris13.fr

Olivier Lopez

Wei Zhang

Giorgio Santarelli

Yann Le Coq

Frédéric Auguste

Benoît Darquié

Christian Chardonnet

Anne Amy-Klein

Mid-IR frequency measurement using an optical frequency comb and a long-distance remote frequency reference

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

LNE-SYRTE

Ultrahigh resolution molecular spectroscopy is an alternative tool to atoms to perform tests of fundamental physics, such as the search for the non conservation of parity [START_REF] Darquié | Progress toward a first observation of parity violation in chiral molecules by high-resolution laser spectroscopy[END_REF], the temporal stability of the electron-to-proton mass ratio [START_REF] Shelkovnikov | Stability of the Proton-to-Electron mass ratio[END_REF], or the electron electric dipole moment [START_REF] Hudson | Improved measurement of the shape of the electron[END_REF]. Many of these tests rely on the availability of ultrastable and accurate laser sources emitting in the mid-IR where molecules exhibit rovibrational transitions. It is thus very challenging to develop a frequency stabilization or measurement scheme in the mid-IR which does not depend on quite rare frequency secondary references. Such high sensitivity spectroscopy experiments have been performed using an optical frequency comb. In that case, the IR frequency is compared to a high-harmonic of the comb repetition rate with a sum or difference frequency generation in a non linear crystal and a Ti:Sa optical frequency comb. It has been demonstrated at 10 µm with a CO 2 laser [START_REF] Amy-Klein | Absolute frequency measurement of an SF6 two-photon line using a femtosecond optical comb and sumfrequency generation[END_REF], at 3.4 µm using a HeNe laser [START_REF] Foreman | Demonstratin of a HeNe/CH4-based molecular clock[END_REF], and recently extended to quantum cascade lasers at 5 µm [START_REF] Borri | Direct link of a mid-infrared QCL to a frequency comb by optical injection[END_REF], and 3.4 µm using a fiber comb [START_REF] Gubin | Femtosecond fiber laser based methane optical clock[END_REF].

II. EXPERIMENTAL SET-UP

In our case, we have built a frequency chain with two fiber femtosecond lasers, a remote optical reference and an optical link, which enables to measure the absolute frequency of a CO 2 laser, emitting in the 28-31 THz frequency range, against an optical reference from the LNE-SYRTE (Fig. 1). The CO 2 laser, located at LPL, is stabilized onto a molecular absorption line.

The optical reference is a fiber laser, emitting at 1.55 µm, locked to a very high finesse cavity with the Pound-Drever-Hall method. A fiber femtosecond laser, centered at 1.55 µm, is phase-locked to the ultrastable laser. The repetition rate of this femtosecond laser is then compared to the primary standards of LNE-SYRTE (H-maser, cryo-oscillator, caesium fountain). This scheme enables to measure the absolute frequency of the ultrastable laser and therefore compensate the cavity drift, which is critical for some high sensitive spectroscopy experiments, as for example the parity violation experiment [START_REF] Darquié | Progress toward a first observation of parity violation in chiral molecules by high-resolution laser spectroscopy[END_REF]. This optical reference is then sent through a 43-km long optical link to LPL [START_REF] Lopez | Cascaded multiplexed optical link on a telecommunication network for frequency dissemination[END_REF].

At LPL, we aim at measuring ?? the frequency of a CO 2 laser, at around 10 µm, locked to an absorption line of the molecule OsO 4 . A fiber femtosecond laser at 1.55 µm is used to transfer the optical reference stability and accuracy to the IR laser. The modes frequencies can be written as:

CEO rep p f f p + = ν (1)
where p is an integer near 780 000. The optical reference is used to lock the repetition rate of the comb. After removing the offset frequency f CEO , we have:

) ( 1 ref ref rep N f Δ - = ν (2)
with f rep the repetition rate, N an integer near 780 000, ν ref the optical reference frequency, and Δ ref the beatnote between the ultrastable laser and the femtosecond laser. Then, lowfrequency modes are generated in a non linear fiber in order to obtain a comb centered around 1.85 µm. This comb is mixed with the CO 2 laser in a non linear crystal leading, via sum frequency generation, to a shifted comb of frequencies:

2 2 CO CEO rep CO q f f q ν ν ν + + = + (3)
with ν q the frequency of one mode of the comb, ν CO2 the CO 2 laser frequency, q an integer around 660 000. The beatnote Δ CO 2 between this shifted comb and the comb directly issued from the femtosecond laser, with modes p ν gives:

rep CO CEO rep CO CEO rep CO f q p f f p f f q ) ( ) ( ) ( 2 2 2 - - = + - + + = Δ ν ν (4)
In combination with (2), we obtain:

2 2 ) ( CO ref ref CO N q p Δ - Δ - - = ν ν (5)
with p-q is near 120 000. Thus the CO 2 laser frequency is compared to a high-harmonic of the comb repetition rate [START_REF] Amy-Klein | Absolute frequency measurement of an SF6 two-photon line using a femtosecond optical comb and sumfrequency generation[END_REF]. Finally, the CO 2 laser frequency can be measured against the frequency standards of LNE-SYRTE.

III. RESULTS AND PERSPECTIVES

With this set-up, we measured the frequency of a CO 2 laser, locked to an OsO 4 saturation line, with a relative stability of 6x10 -14 after 1 s integration time, not limited by the measurement set-up [START_REF] Bernard | CO2 laser stabilization to 0.1-Hz level using external electrooptic modulation[END_REF], even without compensating the optical link used for the experiments (Fig. 2). The measurements duration is however limited by the chemical stability of OsO 4 molecules, which degrade after a few hours. We have also measured the absolute frequency of the CO 2 laser, and we obtained a value consistent with the previous measurements [START_REF] Acef | Accurate determination of OsO4 absolute frequency grid at 28/29 THz[END_REF]- [START_REF] Amy-Klein | Absolute frequency measurement in the 28-THz spectral region with a femtosecond laser comb and a long-distance optical link to a primary standard[END_REF].

We are now progressing towards the stabilization of the CO 2 laser via the frequency comb. Because molecules degrade over a few hours, the CO 2 laser frequency drifts. First results of stabilization show that this drift has been removed. The next step is to lock the CO 2 laser directly to the optical reference,without prestabilization on a molecular absorption line. This is necessary for the parity non conservation experiment since left-and right-handed molecules lines are alternately measured in a supersonic beam [START_REF] Darquié | Progress toward a first observation of parity violation in chiral molecules by high-resolution laser spectroscopy[END_REF].

In the future, we will extend this stabilization/measurement scheme to quantum cascade lasers. Among several constraints, the ideal molecule for the parity violation experiment must have a rovibrational absorption line in the spectral range of an ultrastable laser. Till now, we were limited to the CO 2 laser, which emits in spectral windows about 1 GHz wide every 30-50 GHz, between 9 and 11 µm. With a quantum cascade laser, we can access a continuous spectral range of 200 GHz, at any wavelength between 2 and 20 µm. This would simplify the optical scheme of any ultrahigh sensitivity IR spectroscopy measurement.

Figure 1 .

 1 Figure 1. Mid-IR frequency measurement principle

Figure 2 .

 2 Figure 2. Allan deviation of the CO2 laser frequency locked to an absorption line of OsO4. In red and green, this is respectively the P(16) CO2 laser line near the 190 OsO4 P(55) line and the R(20) CO2 laser line near the 192 OsO4 R(67) line. In blue and magenta, typical stability of optical link used for the IR frequency measurement experiments. Dashed line : combination of the ultrastable laser and primary standards of LNE-SYRTE stabilities.

The authors acknowledge funding support from ANR, LNE, CNRS and Université Paris 13.


Gaëtan Yvert

Marcela Perrone-Bertolotti

Monica Baciu

PhD Olivier David
email: odavid@ujf-grenoble.fr

Dynamic causal modeling of spatiotemporal integration of phonological and semantic processes: An electroencephalographic study

Integration of phonological and lexico-semantic processes is essential for visual word recognition. Here we used dynamic causal modeling of event-related potentials, combined with group source reconstruction, to estimate how those processes translate into contextdependent modulation of effective connectivity within the temporal-frontal language network. Fifteen healthy human subjects performed a phoneme detection task in pseudo-words and a semantic categorization task in words. Cortical current densities revealed the sequential activation of temporal regions, from the occipital-temporal junction towards the anterior temporal lobe, before reaching the inferior frontal gyrus. A difference of activation between phonology and semantics was identified in the anterior temporal lobe, within the 240-300 ms peristimulus time-window. Dynamic causal modeling indicated this increase of activation of the anterior temporal lobe in the semantic condition as a consequence of an increase of forward connectivity from the posterior inferior temporal lobe to the anterior temporal lobe.

In addition, fast activation of the inferior frontal region, that allowed a feedback control of frontal regions on the superior temporal and posterior inferior temporal cortices, was found to be likely. Our results precisely describe spatio-temporal network mechanisms occurring during integration of phonological and semantic processes. In particular, they support the hypothesis of multiple pathways within the temporal lobe for language processing, where frontal regions would exert a top-down control on temporal regions in the recruitment of the anterior temporal lobe for semantic processing.

Introduction

Performing a language task with visually presented stimuli involves completion of several operations, such as visuo-orthographic and grapheme-to-phoneme decoding, phonological and semantic analyses, with lexical retrieval and word recognition [START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF].

Identifying brain network mechanisms associated to those different components of language, e.g. spatio-temporal processes as revealed by chronometry of brain electrical responses or by functional neuroimaging, is important.

Reading written words generates first a visuo-orthographic analysis that mainly involves occipito-temporal cortical areas (Brodmann areas,BA,18,19,37) [START_REF] Greenblatt | Subangular alexia without agraphia or hemianopsia[END_REF][START_REF] Warrington | Word-form dyslexia[END_REF][START_REF] Damasio | The anatomic basis of pure alexia[END_REF][START_REF] Henderson | Anatomy of posterior pathways in reading: a reassessment[END_REF][START_REF] Cohen | The visual word form area: spatial and temporal characterization of an initial stage of reading in normal subjects and posterior split-brain patients[END_REF][START_REF] Leff | The functional anatomy of single-word reading in patients with hemianopic and pure alexia[END_REF][START_REF] Dehaene | The visual word form area: a prelexical representation of visual words in the fusiform gyrus[END_REF][START_REF] Gaillard | Direct intracranial, FMRI, and lesion evidence for the causal role of left inferotemporal cortex in reading[END_REF]. Associated electroencephalographic (EEG) component is a negative wave above occipito-temporal regions at 150 ms after stimulus presentation (N150) [START_REF] Gros | Automatic grapheme processing in the left occipitotemporal cortex[END_REF][START_REF] Cohen | Specialization within the ventral stream: the case for the visual word form area[END_REF][START_REF] Brem | Neurophysiological signs of rapidly emerging visual expertise for symbol strings[END_REF].

Phonological analysis includes grapheme-to-phoneme conversion [START_REF] Berndt | Functional pronunciation units in English words[END_REF], decoding and retrieval of sound forms of the written words [START_REF] Herbster | Regional cerebral blood flow during word and nonword reading[END_REF][START_REF] Fiez | Effects of lexicality, frequency, and spelling-to-sound consistency on the functional anatomy of reading[END_REF][START_REF] Walter | Grapheme-phoneme processing of single words[END_REF][START_REF] Alario | The role of phonological and orthographic information in lexical selection[END_REF][START_REF] Simon | Orthographic transparency and graphemephoneme conversion: An ERP study in Arabic and French readers[END_REF]. Manipulating phonemes and phonological forms of the words involves the superior temporal gyrus (BA 22) [START_REF] Pugh | Cerebral organization of component processes in reading[END_REF][START_REF] Vandenberghe | Functional anatomy of a common semantic system for words and pictures[END_REF][START_REF] Lurito | Comparison of rhyming and word generation with FMRI[END_REF][START_REF] Billingsley | Functional MRI of phonological and semantic processing in temporal lobe epilepsy[END_REF][START_REF] Vigneau | Meta-analyzing left hemisphere language areas: Phonology, semantics, and sentence processing[END_REF] and the inferior frontal gyrus (pars triangularis and opercularis; BA 44, 45) [START_REF] Zatorre | Lateralization of phonetic and pitch discrimination in speech processing[END_REF][START_REF] Dapretto | Form and content: dissociating syntax and semantics in sentence comprehension[END_REF][START_REF] Mcdermott | A procedure for identifying regions preferentially activated by attention to semantic and phonological relations using functional magnetic resonance imaging[END_REF][START_REF] Seghier | Variability of fMRI activation during a phonological and semantic language task in healthy subjects[END_REF]. Semantic processing requiring lexical retrieval of semantic representations is supposed to take place within inferior frontal (BA 47), middle and inferior temporal regions (BA 21, 37) [START_REF] Billingsley | Functional MRI of phonological and semantic processing in temporal lobe epilepsy[END_REF][START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF][START_REF] Gitelman | Language network specializations: An analysis with parallel task designs and functional magnetic resonance imaging[END_REF]. Semantic processing also involves the anterior temporal lobe [START_REF] Fujimaki | Early neural activation for lexico-semantic access in the left anterior temporal area analyzed by an fMRI-assisted MEG multidipole method[END_REF][START_REF] Binney | The Ventral and Inferolateral Aspects of the Anterior Temporal Lobe Are Crucial in Semantic Memory: Evidence from a Novel Direct Comparison of Distortion-Corrected fMRI, rTMS, and Semantic Dementia[END_REF][START_REF] Visser | Semantic processing in the anterior temporal lobes: a meta-analysis of the functional neuroimaging literature[END_REF], which might convert modality specific information into amodal representations in order to extract all semantic attributes [START_REF] Marinkovic | Spatiotemporal dynamics of modality-specific and supramodal word processing[END_REF][START_REF] Patterson | Where do you know what you know? The representation of semantic knowledge in the human brain[END_REF][START_REF] Lau | A cortical network for semantics: (de)constructing the N400[END_REF].

In a meta-analysis [START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF], arguments for a dual route model of reading were provided. In short, this model assumes that grapho-phonological decoding and phonological analysis on the one hand, and lexico-semantic processing on the other hand, are segregated into two parallel pathways. In a recent functional magnetic resonance imaging (fMRI) study using dynamic causal modeling (DCM) [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF], three potential processing streams from occipital to temporal cortices were described: a ventral lexico-semantic route, a dorsal phonological route and an intermediate route encompassing partially both ventral and dorsal route.

Here, we describe spatio-temporal organization of phonological and semantic processing using EEG recordings performed in 15 healthy participants during two separate language tasks involving phonological or semantic processing of visual stimuli. First, we proceed to group-level source reconstruction [START_REF] Litvak | Electromagnetic source reconstruction for group studies[END_REF] to describe the sequence of activation of brain regions involved. Second, we specifically test the face validity of the dual route-model [START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF] using DCM of event-related potentials (ERPs) [START_REF] David | Dynamic causal modeling of evoked responses in EEG and MEG[END_REF][START_REF] Kiebel | Dynamic causal modelling of evoked responses in EEG/MEG with lead field parameterization[END_REF].

Materials and Methods

Participants

Fifteen French volunteers (4 females, 11 males; 23-26 yo, mean 24.5 yo) underwent EEG recording. All participants were right-handed as assessed by the Edinburgh Handedness Inventory [START_REF] Oldfield | The assessment and analysis of handedness: the Edinburgh inventory[END_REF] and had normal or corrected-to-normal vision and no previous neurological episodes or history of substances abuse. They gave informed written consent to the study, which was approved by the Ethics Committee of Grenoble University Hospital (N° ID RCB: 2010-A01063-36).

Stimuli and Tasks

Participants were asked to perform two language tasks during one EEG session: phoneme detection in pseudo-words [START_REF] Perrone-Bertolotti | Functional MRI evidence for modulation of cerebral activity by grapheme-to-phoneme conversion in French, and by the variable of gender[END_REF] and semantic categorization [START_REF] Cousin | Functional MRI approach for assessing hemispheric predominance of regions activated by a phonological and a semantic task[END_REF] in words presented in visual modality. Both tasks were performed separately in two distinct sessions, which were counterbalanced across subjects, in order to avoid needing to distinguish words and pseudowords that would have involved confounding processes at the lexical level. Across both tasks, we controlled the number of letters and syllables, as well as the bigram and trigram frequencies according to Lexique 2 database.

Phoneme detection

Stimuli were 192 legal French pseudo-words (5-7 letters length) built by exchanging two or three letters in French concrete nouns. Half of the items were considered "target" and contained the target phoneme to be detected (/o/). The grapheme-to-phoneme conversion (GPC) for the phoneme /o/ had two levels of difficulty. For half of the target items, the grapheme representation was simple (Target S) as phoneme /o/ corresponded to grapheme ''o" (e.g., obalet, clotir, damulo). For the other half, the grapheme representation of phoneme /o/ was complex (Target C) as the sound /o/ corresponded to graphemes "au" (e.g., autril, phaumi, amilau). Furthermore, the phoneme /o/ was presented an equivalent number of times at the beginning (e.g., autril or obalet), in the middle (e.g., phaumi or clotir) and at the end (e.g., amilau or damulo) of the pseudo-words. Participants were asked to judge whether the pseudo-words contained the sound /o/ or not. They were instructed to internally pronounce the item, without articulation and verbalization.

Semantic categorization

Stimuli were 192 words (5-7 letters length). Half of the items were "target" and described "living" entities (plants and animals). The other half represented "non living" entities. They were selected from Lexique 2 database for French [START_REF] New | Lexique 2: a new French lexical database[END_REF]. Lexical frequencies of words and the number of letters (between 5 and 7) were controlled between both conditions.

Participants were asked to judge whether the word described a living or non-living entity.

They were instructed to read words without articulation and verbalization.

Stimuli were generated using E-Prime software running on a PC (E-prime Psychology Software Tools Inc., Pittsburgh, USA). Stimuli were written in white "Courier New", font size 24, centered on a middle of a black screen and presented for 1000 ms. A fixation cross was presented between stimuli. Participants performed manual responses by pressing two keys on a keyboard with the index and the middle finger of the right hand. Each trial lasted 2.5 sec and consisted of 500 ms fixation cross, 1000 ms pseudo-word/word presentation and 1000 ms fixation cross for providing the response. The reaction time (RT) and the accuracy (% correct responses, % CR) were recorded for each participant in order to evaluate the performance of task execution. Each task lasted approximately 8 min.

EEG acquisition

After positioning the electrodes on the scalp, participants were seated in a comfortable chair in front of the projection screen (70 cm distance). During sessions and between trials, they were asked to fixate a cross in the middle of the screen and to refrain from moving and blinking. EEG data were recorded using 96 Ag-AgCl active electrodes (BrainCap™, Brain Products GmBH, Gilching, Germany) positioned according to the standardized 10-20 system. EEG signals were recorded continuously with a sampling rate of 500 Hz and electrode impedances were kept below 5 kΩ. Reference and ground were taken at FCz (midline, central frontal) and AFz (midline, anterior frontal), respectively. Positions of the 96 electrodes were obtained as follows: the three fiducials (nasion, left ear, right ear) and six well-distributed electrodes (Fp1, Fp2, P3, P4, Cz, Oz) were measured manually with callipers [START_REF] Koessler | Spatial localization of EEG electrodes[END_REF]. Those measures were then used to obtain the positions of all electrodes after computation of the affine transformation that best fitted the canonical EEG cap provided in SPM8 onto the head of each subject (Wellcome Trust Centre for Neuroimaging, http://www.fil.ion.ucl.ac.uk/spm/).

EEG processing

EEG data were processed using SPM8 for EEG software (update 4290).

EEG pre-processing

Continuous EEG data were down-sampled to 250 Hz, band-pass filtered between 1 and 30 Hz (Butterworth filter) and converted to average reference. They were then epoched and baseline-8 corrected from -500 ms to 1000 ms, where time 0 is stimulus onset. Trials including eyes movements (>80 mV), muscular artifacts (high frequencies), and poor responses (wrong answer or reaction time of more than a second) were excluded from further analysis (on average over subjects, for phonologic and semantic tasks respectively, 84.43 % ± 9.33 % and 82.67% ± 8.86 % of trials were kept). ERP were obtained by averaging epoched trials for each subject and each task.

Statistical analyses of ERP

Statistical analyses of ERP at the group level were performed according to the SPM procedure (Kilner and Friston, 2010). First, ERP amplitudes for each subject and each condition were transformed into scalp maps of dimension 64x64 obtained for the 375 time bins using a twodimensional (2D) linear interpolation (knowing electrode coordinates in the Montreal Neurological Institute space -MNI). Second, scalp-maps of ERP amplitudes were concatenated over time to produce a three-dimensional (3D) volume (2D space x time). Third, the scalp-time 3D volume was smoothed using a low-pass kernel ([9 mm x 9 mm x 20 ms] full width at half magnitude -FWHM) to obtain data that could be approximated to a Gaussian random field in the absence of effect. Group statistics on this scalp-time volume were finally obtained using a one sample t-test to detect common responses to both conditions and a paired t-test to detect differences between conditions. Statistical analyses were performed from -500 ms up to 1000 ms and were corrected for multiple comparisons by controlling the family wise error (FWE) rate across space and time.

Cortical source imaging

Task-related cortical current source densities were obtained from ERP using an inversion procedure. A canonical mesh (cortex, skull, skin) was used to ensure that cortical activity is reconstructed in the same source space over subjects [START_REF] Mattout | Canonical Source Reconstruction for MEG[END_REF]. The meshes were warped to each subject's anatomy using an inverse spatial normalization defined with the three fiducials and the electrode coordinates. The canonical cortical mesh was composed of 8196 vertices. A boundary element method (BEM) head model composed of three layers (inner skull, outer skull and scalp with respective conductivities of: 0.33 S.m -1 , 0.0041 S.m -1 and 0.33 S.m -1 ) was then used to calculate the forward model that links cortical activity to scalp measurements [START_REF] Phillips | Forwards models for EEG[END_REF]. Cortical current source densities were obtained using the group inversion procedure described in [START_REF] Litvak | Electromagnetic source reconstruction for group studies[END_REF]. It is based on a hierarchical model with multiple sparse priors for bioelectric source imaging [START_REF] Friston | Multiple sparse priors for the M/EEG inverse problem[END_REF] with a particular hyperprior that imposes consistency on the source reconstruction over conditions and subjects. In other words, source reconstruction in each individual analysis was confined to the same locations and only amplitudes of cortical current densities were allowed to be modulated by the task. The temporal window of inversion was set between 200 ms prestimulus and 600 ms post-stimulus. Two inverse solutions at the group level were performed:

(i) pooling together both experimental conditions to focus on differences of amplitude of activation; (ii) separating experimental conditions to focus on differences of source location.

The second option was only chosen to optimise the choice of regions of interest (ROI) used in DCM (see below) after finessing the source analysis in the 140-240 ms range, which contains the second ERP component corresponding to early processing of both tasks in occipitaltemporal regions. For every subject and condition, 3D images in the MNI space of the energy of current density was obtained every 20 ms by applying a temporal Gaussian kernel of 20 ms FWHM to cortical time series and an isotropic Gaussian kernel of 9 mm FWHM used for spatial resampling.

Group-level T-statistics of 3D images of cortical energy were obtained using a standard random-effect analysis in SPM8. The image of every post-stimulus time-window was contrasted against a baseline time window located between -100 ms and -80 ms peristimulus time. Small volume correction for multiple comparisons was applied into a 10 mm radius sphere centered on the cluster of activation (p<0.005 uncorrected). Final statistical threshold was set at p<0.05, FWE corrected for multiple comparisons at the cluster level.

Dynamic causal modeling

DCM for ERP explains ERP as a resultant of stimulus-related changes in the activity of neuronal populations [START_REF] David | Dynamic causal modeling of evoked responses in EEG and MEG[END_REF][START_REF] Kiebel | Dynamic causal modelling of evoked responses in EEG/MEG with lead field parameterization[END_REF]. Each cortical source is represented by three interconnected populations (pyramidal cells, excitatory interneurons and inhibitory interneurons) representing different cortical layers. Different sources are connected by long-range connections conforming to physiological connectivity rules [START_REF] Felleman | Distributed hierarchical processing in the primate cerebral cortex[END_REF]. Differences in interregional effective connectivity afford changes in the shape of ERP between conditions. The models and their constitutive sources are first specified according to functional hypotheses of interaction across brain structures of interest. Second, model parameters are inferred from the ERP and the evidences of models are computed using Bayesian procedures [START_REF] Penny | Comparing dynamic causal models[END_REF]. Finally, Bayesian model selection (BMS) is applied to define the most plausible model [START_REF] Stephan | Bayesian model selection for group studies[END_REF], or family of models [START_REF] Penny | Comparing families of dynamic causal models[END_REF]. It computes the expected probability of obtaining a given model for any randomly selected subject in the group, and the exceedance probability, which is the belief that a particular model is more likely than any other model (of all models tested), given the group data. The BMS procedure takes into account not only the goodness of fit of the data, but also the complexity of the models (i.e. number of free parameters) in order to avoid overfitting the data. Averaged parameters of relatively similar, but different, models can be obtained or each subject using Bayesian model averaging (BMA), which averages the parameters of the different models weighted by their posterior probability [START_REF] Trujillo-Barreto | Bayesian model averaging in EEG/MEG imaging[END_REF].

Inference at the group level on connectivity parameters obtained after BMA of selected models can then be performed to identify significant changes of connectivity between conditions.

Dynamic causal models (DCMs) were constructed with five ROIs identified from our source localization results (Table 1, see also Figures 3 and4): occipito-temporal junction (OT), posterior part of inferior temporal (IT), posterior part of superior temporal (STS), anterior temporal (AT) and inferior frontal gyrus (IFG). Each ROI was modeled with an equivalent current dipole (ECD) [START_REF] Kiebel | Dynamic causal modelling of evoked responses in EEG/MEG with lead field parameterization[END_REF], which was positioned symmetrically in each hemisphere on the maximum of activation clusters obtained after group source reconstruction (see Table 1 for MNI coordinates). Note that we did not explicitly include early relays within the visual thalamus and primary visual cortex to keep the models as simple as possible. This early activity is implicitly modeled in the shape of the extrinsic input to OT, which is estimated from the data, assuming identical early visual processing between conditions. All models included a ventral pathway (OT / IT / AT), a dorsal pathway (OT / STS / AT) and reciprocal connections between IT and STS [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF]. They were symmetric across hemispheres, including lateral non task-modulated connections between homologous regions. Sixteen models were designed in order to test three hypotheses (Figure 1): (i) the level in the cortical hierarchy of connectivity modulation within the dorsal and ventral routes (just before or just after STS / IT); (ii) the presence of a fast pathway towards IFG that could facilitate the top-down influence of IFG on lower regions [START_REF] Cornelissen | Activation of the Left Inferior Frontal Gyrus in the First 200 ms of Reading: Evidence from Magnetoencephalography (MEG)[END_REF]; (iii) the feedback pathway of IFG on STS or IT or both. These sixteen models were computed with and without modulation of intrinsic connections, i.e. between neuronal populations within a given ROI, resulting in the study of 32 models. The evidence of all models was assessed on ERP time series taken between 0 and 450 ms peri-stimulus time. DCM options when adjusting data were chosen as follows: eight modes for data reduction; one DCT component to remove slow drift; Hanning windowing to ensure convergence and limit the effect of late components (after 350 ms); no constraint on the orientation symmetry of ECD to favor differences of activation between left and right hemispheres. Models were compared at the group level using random effect [START_REF] Stephan | Bayesian model selection for group studies[END_REF]. Families of models [START_REF] Penny | Comparing families of dynamic causal models[END_REF] were also compared to test four hypotheses (Figure 1

Results

Behavior

Statistical analyses of accuracy (% of correct responses, %CR) and reaction times (RT, ms) revealed that phonological (M = 84.43 %; SD = 9.33 %) and semantic (M = 82.67 %; SD = 8.86 %) tasks were performed correctly. No significant difference (p = 0.8) was obtained between tasks in terms of accuracy. However, the speed of execution was significantly higher (T(14) = 3.7; p<0.01) for the semantic task (M = 673.34 ms; SD = 50.55 ms) than for the phonological task (M = 720.25 ms; SD = 32.74 ms). Because basic linguistic features (sublexical level), i.e. number of letters, phonemes and frequencies of bigrams and trigrams, were equivalent for both tasks, the difference of execution speed is more likely to have been related to high level processes, i.e. lexical or semantic.

ERP

Root mean square (RMS) of the ERP amplitude averaged over subjects and electrodes showed for each condition three main components occurring at 112 ms (P100), 196 ms (N170/N200) and 332 ms (P350) with respect to stimulus onset (Figures 2A and2B).

Cortical source imaging

In order to find the common network for both phonological and semantic tasks, current source densities were first estimated using ERP obtained for both conditions, from 100 to 350 ms (Figure 3A). The P100-140 observed on occipital electrodes was estimated to originate from the right inferior and middle occipital gyri. The N150 was located in the left inferior occipitaltemporal cortex. From 160 to 220 ms, the positive ERP was located across superior, middle and inferior temporal gyri, bilaterally but more predominant to the left. From 220 ms to 260 ms, the activations spread across the temporal lobe from the inferior gyrus to the anterior pole, bilaterally. After 260 ms, the activations spread from the temporal pole towards the inferior frontal gyrus bilaterally. The difference of cortical activation between both conditions was marginally significant (p<0.001 uncorrected) between 240 and 300 ms only. It was located within the anterior part of the temporal lobe, predominantly to the left, and corresponded to higher responses of this region during the semantic task (Figure 3B).

To clarify whether feedforward information transfer through the temporal lobe would have involved, not only different activation strength, but also different regions between semantic and phonological processing, ERP were also group-inverted using the same methodology, but separately for each condition, between 140 and 240 ms (

Dynamic Causal Modeling

Bayesian model selection at the model level indicated Model 14 as being the most plausible one (Figure 5A). This model included a fast connection from OT to IFG, a feedback from IFG to STS and a late modulation of temporal connectivity. Interestingly, the selection of this model could not be attributed solely to its complexity, because complexity of Model 8 and Model 9 or 15 was equivalent and most of the models with intrinsic modulation (Models 1732) were more complex than Model 14. Although comparison of model families clearly indicated a modulation between higher temporal regions rather than between lower ones (Figure 5B) and a direct connection from OT to IFG (Figure 5C), in accordance with the winning model, comparison of families depending on the target of the feedback was not conclusive except that a feedback was strongly needed, either on STS, IT, or STS/IT (Figure 5D). Finally, models with or without modulation of intrinsic connectivity were equally likely (Figure 5E).

According to these results, we considered the family of 6 models (Models 14-15-16 and their homologues with intrinsic modulation: Models 30-31-32) that comprised the most likely features (anterior level of modulation, connection between OT and IFG, feedback from IFG)

and applied the BMA procedure to obtain connectivity parameters of that family for both conditions. A paired t-test on the log connectivity parameters revealed a significant increase (p<0.05, Bonferroni corrected) of the connectivity from IT to AT in the semantic task compared to the phonological task (Figure 6). In addition, we found an increase of connectivity between IT and STS in the semantic condition compared to the phonological one. It was not significant (p=0.17, Bonferroni corrected) when testing log connectivity parameters estimated from BMA, but reached significance (p<0.05, Bonferroni corrected) when testing log connectivity parameters averaged across models using the standard averaging procedure, i.e. without weighting by the posterior probability of models.

On average, cortical time series did not differ a lot between conditions. Only a marginal effect The aim of this EEG study was to identify the spatio-temporal organization of brain response to word recognition by using phonological and semantic tasks. We based our work on the assumptions of the dual route model of word reading [START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF] which considers two specific pathways, one for phonological and another one for semantic processing. It is in line with a recent DCM study in fMRI showing multiple routes from occipital to temporal cortex during reading [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF]. Methodologically, we measured the neural responses using group source reconstruction [START_REF] Litvak | Electromagnetic source reconstruction for group studies[END_REF]) and DCM analyses [START_REF] David | Dynamic causal modeling of evoked responses in EEG and MEG[END_REF].

Our main results suggest different pathways for word recognition and processing. One pathway, that links the occipital-temporal junction to the anterior temporal cortex via ventral temporal regions, would be predominantly required for the semantic task. Another pathway, also linking the occipital-temporal junction to the anterior temporal cortex, but going through dorsal temporal regions, would be required more for lexico-phonological processing. In addition, DCM indicated: (i) fast activation of inferior frontal cortex suggesting top-down effects on the posterior part of the temporal lobe; (ii) increase of forward connectivity from inferior posterior to anterior and superior temporal regions in the semantic task compared to the phonological task.

Spatio-temporal integration of phonological and semantic processes

Group source reconstruction showed activation of a predominant left cerebral network classically involved in language processing and including left OT, left posterior IT, left posterior STS, bilateral AT and bilateral IFG.

Early activation of OT suggests sub-lexical orthographic processing in that region (Maurer et al., 2005;[START_REF] Dien | The neurocognitive basis of reading single words as seen through early latency ERPs: a model of converging pathways[END_REF]. OT has been reported to be sensitive to letter string analysis [START_REF] Cohen | The visual word form area: spatial and temporal characterization of an initial stage of reading in normal subjects and posterior split-brain patients[END_REF][START_REF] Tarkiainen | Dynamics of visual feature analysis and object-level processing in face versus letter-string perception[END_REF] and could be involved in the segmentation of word and pseudo-word into bigram and trigram entities (graphemes) [START_REF] Dehaene | The visual word form area: a prelexical representation of visual words in the fusiform gyrus[END_REF][START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF].

Classically, increased complexity of stimulus processing is associated with posterior-toanterior brain gradient according to the change in stimulus representation, from unimodal to multimodal and amodal [START_REF] Nobre | Word recognition in the human inferior temporal lobe[END_REF][START_REF] Buchel | A multimodal language region in the ventral visual pathway[END_REF][START_REF] Binder | Where Is the Semantic System? A Critical Review and Meta-Analysis of 120 Functional Neuroimaging Studies[END_REF][START_REF] Dien | The neurocognitive basis of reading single words as seen through early latency ERPs: a model of converging pathways[END_REF]. In this framework, during phonological task, rapid activation (160-200ms) of the posterior STS involved in multimodal representation [START_REF] Beauchamp | Integration of auditory and visual information about objects in superior temporal sulcus[END_REF][START_REF] Reale | Auditory-visual processing represented in the human superior temporal gyrus[END_REF], may account for direct grapheme-to-phoneme conversion [START_REF] Simos | Brain mechanisms for reading: the role of the superior temporal gyrus in word and pseudoword naming[END_REF]. During semantic task, rapid activation (160-220 ms) of the posterior part of inferior temporal gyrus, also involved in multimodal representation [START_REF] Price | The myth of the visual word form area[END_REF][START_REF] Dien | The neurocognitive basis of reading single words as seen through early latency ERPs: a model of converging pathways[END_REF][START_REF] Kherif | Automatic Top-Down Processing Explains Common Left Occipito-Temporal Responses to Visual Words and Objects[END_REF], may be related to graphemic to lexical matching within mental lexicon [START_REF] Martin | The Representation of Object Concepts in the Brain[END_REF][START_REF] Lau | A cortical network for semantics: (de)constructing the N400[END_REF]. The rapid segregation of phonological and lexical pathways is in line with results reported in [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF]. From 180 to 220 ms, co-activation of posteroinferior and superior temporal regions would be related to lexico-phonological representations of words [START_REF] Wise | Separate neural subsystems within 'Wernicke's area[END_REF][START_REF] Hickok | The cortical organization of speech processing[END_REF][START_REF] Gagnepain | Spoken Word Memory Traces within the Human Auditory Cortex Revealed by Repetition Priming and Functional Magnetic Resonance Imaging[END_REF].

In the present study, the anterior part of the inferior temporal lobe was subsequently activated for both tasks, from 240 to 280 ms, in agreement with similar results provided by a MEG study [START_REF] Fujimaki | Early neural activation for lexico-semantic access in the left anterior temporal area analyzed by an fMRI-assisted MEG multidipole method[END_REF]. However, our results revealed longer involvement of the activation of the anterior part of inferior temporal lobe for semantic, compared to phonological task. This result confirms the crucial role of this region for the semantic access [START_REF] Mummery | A voxel-based morphometry study of semantic dementia: relationship between temporal lobe atrophy and semantic memory[END_REF][START_REF] Marinkovic | Spatiotemporal dynamics of modality-specific and supramodal word processing[END_REF][START_REF] Bright | Unitary vs multiple semantics: PET studies of word and picture processing[END_REF][START_REF] Gorno-Tempini | Cognition and anatomy in three variants of primary progressive aphasia[END_REF][START_REF] Noppeney | Temporal lobe lesions and semantic impairment: a comparison of herpes simplex virus encephalitis and semantic dementia[END_REF][START_REF] Fujimaki | Early neural activation for lexico-semantic access in the left anterior temporal area analyzed by an fMRI-assisted MEG multidipole method[END_REF][START_REF] Ralph | Conceptual knowledge is underpinned by the temporal pole bilaterally: convergent evidence from rTMS[END_REF][START_REF] Binney | The Ventral and Inferolateral Aspects of the Anterior Temporal Lobe Are Crucial in Semantic Memory: Evidence from a Novel Direct Comparison of Distortion-Corrected fMRI, rTMS, and Semantic Dementia[END_REF][START_REF] Visser | Semantic processing in the anterior temporal lobes: a meta-analysis of the functional neuroimaging literature[END_REF]. This region could serve as a hub within the widespread semantic system. Its role would be to convert lexical into amodal representations [START_REF] Marinkovic | Spatiotemporal dynamics of modality-specific and supramodal word processing[END_REF][START_REF] Patterson | Where do you know what you know? The representation of semantic knowledge in the human brain[END_REF][START_REF] Lau | A cortical network for semantics: (de)constructing the N400[END_REF][START_REF] Pulvermuller | The word processing deficit in semantic dementia: all categories are equal, but some categories are more equal than others[END_REF].

Finally, the left IFG was robustly and equally activated for both tasks from 280 to 340 ms.

The activity in this region may suggest top-down control of semantic memory in order to guide retrieval and selection of conceptual information stored within posterior temporal cortex [START_REF] Thompson-Schill | Verb generation in patients with focal frontal lesions: a neuropsychological test of neuroimaging findings[END_REF][START_REF] Badre | Dissociable controlled retrieval and generalized selection mechanisms in ventrolateral prefrontal cortex[END_REF].

Dynamic causal modeling of the written language network

DCM was used here to further address effective connectivity underlying the observed sequence of activation. In particular, we concentrated on the multiple routes model [START_REF] Jobard | Evaluation of the dual route theory of reading: a metanalysis of 35 neuroimaging studies[END_REF][START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF], using the ROIs discussed just above. A ventral pathway (OT / IT / AT) and a dorsal pathway (OT / STS / AT) were considered, each one supposed to segregate semantic and phonological processing, respectively. Those regions also embedded an intermediate lexico-phonologic pathway (OT / IT / STS / AT). Specifically, 32 DCM models were tested to address specific hypotheses through Bayesian model selection: (i) the level in the cortical hierarchy of connectivity modulation within the dorsal and ventral routes;

(ii) the presence of a fast pathway towards IFG that could facilitate the top-down influence of IFG on lower regions; (iii) the feedback pathway of IFG on STS or IT or both. These models recapitulate to some extent the pathways from basic visual to higher-order temporal language areas proposed in [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF]: (1) OT / IT / AT; (2) OT / STS / AT; (3) OT / IT / STS / AT.

We detected significant modulation of effective connectivity on the forward connection between IT and AT, which was increased during the semantic condition. These regions are connected anatomically via the longitudinal fasciculus [START_REF] Crosby | Correlative Anatomy of the Nervous System[END_REF][START_REF] Catani | Occipito-temporal connections in the human brain[END_REF] and are functionally included in the pathway (1) OT / IT / AT proposed by [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF]. This result is in agreement with models suggesting that written words can access semantics independently of phonology [START_REF] Seidenberg | A distributed, developmental model of word recognition and naming[END_REF][START_REF] Harm | Computing the meanings of words in reading: cooperative division of labor between visual and phonological processes[END_REF]) by direct graphemic-lexical matching within IT [START_REF] Lau | A cortical network for semantics: (de)constructing the N400[END_REF] and amodal representation within AT [START_REF] Marinkovic | Spatiotemporal dynamics of modality-specific and supramodal word processing[END_REF][START_REF] Patterson | Where do you know what you know? The representation of semantic knowledge in the human brain[END_REF].

Furthermore, increased connectivity between IT and STS during semantic condition may account for the proposed pathway (3) OT / IT / STS / AT of [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF], which possibly reflects the integration between lexical and phonological representations of the word [START_REF] Wise | Separate neural subsystems within 'Wernicke's area[END_REF][START_REF] Hickok | The cortical organization of speech processing[END_REF][START_REF] Gagnepain | Spoken Word Memory Traces within the Human Auditory Cortex Revealed by Repetition Priming and Functional Magnetic Resonance Imaging[END_REF]. The increase of connectivity between IT and STS was however marginal, thus reflecting a certain degree of inter-subject variability of word processing within the temporal lobe. The pathway (2) OT / STS / AT proposed by [START_REF] Richardson | Multiple Routes from Occipital to Temporal Cortices during Reading[END_REF], which may be purely phonological was not revealed by our analysis. This discrepancy could be due to the fact that our tasks could have equally activated this pathway, particularly the connection between STS and AT.

The comparison between models revealed robust and direct connection from OT to IFG. It allows early activation of frontal cortex during visual language [START_REF] Cornelissen | Activation of the Left Inferior Frontal Gyrus in the First 200 ms of Reading: Evidence from Magnetoencephalography (MEG)[END_REF] and object recognition [START_REF] Bar | A cortical mechanism for triggering top-down facilitation in visual object recognition[END_REF][START_REF] Bar | Top-down facilitation of visual recognition[END_REF] tasks, possibly by the means of direct anatomical connections between early visual regions and frontal cortex [START_REF] Barbas | Connections underlying the synthesis of cognition, memory, and emotion in primate prefrontal cortices[END_REF][START_REF] Kringelbach | The functional neuroanatomy of the human orbitofrontal cortex: evidence from neuroimaging and neuropsychology[END_REF][START_REF] Makris | Segmentation of Subcomponents within the Superior Longitudinal Fascicle in Humans: A Quantitative, In Vivo, DT-MRI Study[END_REF][START_REF] Wakana | Fiber tract-based atlas of human white matter anatomy[END_REF][START_REF] Bürgel | White matter fiber tracts of the human brain: Three-dimensional mapping at microscopic resolution, topography and intersubject variability[END_REF].

Functionally, visual stimuli may activate very early IFG in order to prime the language network for the upcoming cross-modal interaction between visual and language systems for word recognition [START_REF] Cornelissen | Activation of the Left Inferior Frontal Gyrus in the First 200 ms of Reading: Evidence from Magnetoencephalography (MEG)[END_REF]. Additionally, IFG may exert feedback control on regions involved in lower level analysis of the written words [START_REF] Gold | Dissociation of Automatic and Strategic Lexical-Semantics: Functional Magnetic Resonance Imaging Evidence for Differing Roles of Multiple Frontotemporal Regions[END_REF][START_REF] Kherif | Automatic Top-Down Processing Explains Common Left Occipito-Temporal Responses to Visual Words and Objects[END_REF]. However, the precise role and the target region of this feedback remain unclear. The feedback may also reflect attentional processes to control integration of bottom-up processes by lower-level analysis regions depending on the task demand [START_REF] Thompson-Schill | The frontal lobes and the regulation of mental activity[END_REF][START_REF] Gold | Dissociation of Automatic and Strategic Lexical-Semantics: Functional Magnetic Resonance Imaging Evidence for Differing Roles of Multiple Frontotemporal Regions[END_REF][START_REF] Twomey | Top-down modulation of ventral occipito-temporal responses during visual word recognition[END_REF]. The target of the feedback could also depend on the strategy adopted by the subject. Indeed, particularly in the phonological task, subjects could have adopted a phonological strategy, which would be implemented with a feedback on STS, but they could also have used a more orthographic strategy by detecting the target at the grapheme level with a feedback on IT. Further investigations are needed to disentangle this issue. Another limitation of this study regarding IFG is that we were unable to reproduce the late activation of IFG with the tested DCMs. Because we also identified that right parietal activity was also poorly fitted at the scalp level, we can hypothesize that, although tested DCMs were well suited to model fast activation of temporal regions (scope of this study), an additional dorsal pathway could have been required to facilitate the fit of late (>300 ms) frontal activity, e.g. by the means of reciprocal parietal connections with posterior temporal regions and IFG. It would be interesting to clarify this issue in future studies.

Figures & Tables

  /2124/2932); (ii) OT/IFG connectivity (Models 18 and 16 24 vs. Models 916 and 2532); (iii) position of IFG feedback(without feedback Models 1-5-9-13-17-21-25-29 vs. feedback on STS Models 2-6-10-14-18-22-26-30 vs. feedback on IT Models 3-7-11-15-19-23-27-31 vs. feedback on both IT and STS Models 4-8-12-16-20-24-28-32). (iv) presence of intrinsic modulation (Models 116 vs. Models 1732). ********************************* Insert Table 1 & Figure 1 about here ********************************

  Tracking of ERP changes over time at the group level revealed propagation of the activity from occipital towards temporal electrodes and before reaching left frontal regions, for both language tasks: The P100 was detected on the posterior electrodes, predominantly in the right hemisphere in the vicinity of the occipital electrode PO8. The N170/N200 spread from left temporal-occipital electrodes (160 ms) towards more central region (200 ms). The cortical wave of activation finally converged in bilateral temporal regions (240 ms). The maximum of the P350 was detected on left frontal electrodes (320-340 ms). Paired t-tests between task conditions revealed statistical differences only between 240-300 ms on left temporal electrodes [T(14) = 4.8, p <0.001 uncorrected] (Figure2C). Time series averaged over those electrodes showed during this time period a difference in the slope of the ERP between the N200 and P350 components (Figure2D). ********************************** Insert Figure 2 about here ********************************

  Figure 4): Between 140-160 ms, both tasks activated the inferior part of the occipital-temporal cortex; between 160-180 ms, whereas activation for the phonological task moved to the upper part of the posterior temporal cortex, activation for the semantic task remained in the posterior part of the inferior temporal cortex; Between 180-220 ms, activations moved towards more anterior regions, for the phonological task in the upper bank of the temporal cortex and for the semantic task in the posterior part of the temporal cortex. Finally, around 220-240 ms, activation of both tasks converged in the inferior part of the anterior temporal cortex. This result suggested subtle differences in early temporal pathways between the two tasks that were included in the DCMs. ********************************** Insert Figures 3 & 4 about here **********************************

(

  p<0.1 uncorrected) was found at the level of left AT around 240 ms. In addition, a strong left predominance was mostly found in AT. Note also that DCM IFG time series did not show any strong activity after 250 ms, as opposed to cortical current densities. A post-hoc analysis revealed that the fit of the data by DCM was excellent, except after 280 ms where data on few left frontal and right parietal electrodes were not well fitted (results not shown). ********************************** Insert Figures 5 & 6 about here ********************************** 4 Discussion

Figure 1 .Figure 2 :

 12 Figure 1. List of DCMs (M1 to M32). Stimulus extrinsic input is on OT. All connections

Figure 3 .

 3 Figure 3. Source reconstruction averaged over 20 ms time windows from 100 ms to 340 ms.

Figure 4 .

 4 Figure 4. Source reconstruction averaged over 20 ms time windows of the phonological

Figure 5 .

 5 Figure 5. Bayesian model selection between the 32 models (column A) and between families

Figure 6 .

 6 Figure 6. Winning model family composed of a direct connection from OT to IFG, a late

Table 1 :

 1 ROIs used for DCM. MNI coordinates indicate the positions of equivalent current dipoles.

	MNI coordinates

We thank Dr. Vladimir Litvak for his methodological support. This study was funded by Région Rhône-Alpes.