Carine André
email: carine.andre@lpl.univ-aix.fr

Alain Ghio
email: alain.ghio@lpl.univ-aix.fr

Christian Cavé
email: christian.cave@lpl.univ-aix.fr

Bernard Teston
email: bernard.teston…@lpl.univ-aix.fr

PERCEVAL: a Computer-Driven System for Experimentation on Auditory and Visual Perception

Since perception tests are highly time-consuming, there is a need to automate as many operations as possible, such as stimulus generation, procedure control, perception testing, and data analysis. The computer-driven system we are presenting here meets these objectives. To achieve large flexibility, the tests are controlled by scripts. The system's core software resembles that of a lexical-syntactic analyzer, which reads and interprets script files sent to it. The execution sequence (trial) is modified in accordance with the commands and data received. This type of operation provides a great deal of flexibility and supports a wide variety of tests such as auditory-lexical decision making, phoneme monitoring, gating, phonetic categorization, word identification, voice quality, etc. To achieve good performance, we were careful about timing accuracy, which is the greatest problem in computerized perception tests.

INTRODUCTION

Perception tests are time-consuming, not only because several speakers and listeners must be used to compensate for individual variability, but also because measurements must be repeated many times to allow for statistical validation of the results. Hence, it is necessary to computerize in auditory-and visual-perception experiments.

Different systems are available: PsyScope [START_REF] Cohen | PsySCope: a new graphic interactive environment for designing psychology experiments[END_REF], EXPE [START_REF] Pallier | EXPE: an Expandable Programming Language for On-line Psychological Experiments[END_REF], DMDX [START_REF] Forster | DMDX: A Windows Display Program with Millisecond Accuracy[END_REF], Inquisit, etc. Generally, such software supports a wide variety of psychological data-collection methods. Some are easy to use, especially if they provide a graphic interface for experimentation. However, such tools are often designed for specific kinds of experiments derived from experimental psychology, and they are not suitable for other fields, particularly not for the phonetic sciences. In addition, these programs often have to be used with obsolete operating systems or special-purpose hardware (response box/sound board), which are always sources of trouble.

To develop the new version of our experimental package called PERCEVAL (French acronym for PERCeption and EVALuation), we set the following constraints: 1. Ability to simultaneously present visual and audio stimuli. 2. Response recorded through a standard keyboard or a response button box. 3. Response time measured with a temporal resolution of about 1 ms.

4. Maximum flexibility for designing perception tests. 5 . Easy-to-use software with a user-friendly feature control panel. 6. Compatibility with the collective workstation [START_REF] Cavé | A Computer-Driven System for Assessing Speech Quality and Intelligibility[END_REF], on which the system tests up to eight subjects simultaneously. 7. Up-to-date operating system.

SYSTEM CONSIDERATIONS

The basic system has the following components:

•

A standard PC-type computer running under Windows 98/2000 OS. It is equipped with a SVGA graphic screen and a multimedia sound board.

•

A custom-made button box can be added to record the subject's responses. These boxes are built with the printed circuit board (PCB) of a USB joystick or game pad. PCBs are inserted in plastic cases with up to 8 push buttons. These buttons are short stroke, sharp contact, fast action microswitches specially suited to measuring reaction times.

•

For auditory stimulation, we generally use an audio amplifier and Beyer BT100 electrodynamic earphones. They were chosen for their durability, high sensitivity, and response curves [START_REF] Cavé | A Computer-Driven System for Assessing Speech Quality and Intelligibility[END_REF]. The system can be used on a portable PC computer.

TIMING ACCURACY

Timing accuracy is important in perception tests. It can be critical in priming experiments in psychology, for example, where the delay between two views must be precisely controlled. Response time is also important. It can inform us about the difficulty of the task, and it may be the only information when the effect is saturated [START_REF] Mora | Some corpora in (Venezuelian) Spanish for Rapid TTS Evaluation Using a Sentence Verification Task[END_REF].

Why is timing accuracy difficult to control on computer-driven perception tests?

Controlling and measuring time on computer-driven perception tests has been studied for a long time [START_REF] Cohen | PsySCope: a new graphic interactive environment for designing psychology experiments[END_REF][START_REF] Pallier | EXPE: an Expandable Programming Language for On-line Psychological Experiments[END_REF][START_REF] Forster | DMDX: A Windows Display Program with Millisecond Accuracy[END_REF][START_REF] Cavé | A Computer-Driven System for Assessing Speech Quality and Intelligibility[END_REF][START_REF] Macinnes | Millisecond Timing on PCs and Macs[END_REF][START_REF] Howard-Jones | SOAP' -A Speech Output Assessment Package for Controlled Multilingual Evaluation of Synthetic Speech[END_REF]. Programs that run on a single-task operating system (e.g. DOS) have the advantage of being able to perfectly time of the computer-driven test. In fact, in such an environment, the program runs sequentially without "long interruptions" and most of the computer's computational power is dedicated to the test program. On the other hand, with multitask operating systems (e.g. Microsoft Windows), a program may be interrupted for a long time because the computer shares its computational power among all applications running at the same time. In computer-driven perception tests, this phenomenon can be catastrophic if, for example, the dedicated program is inactive when the screen has to be refreshed or when the subject is responding. This can cause timing inaccuracy.

Why is it necessary to develop packages on multitask operating systems like Windows?

Programs developed for DOS generally require dedicated hardware [START_REF] Pallier | EXPE: an Expandable Programming Language for On-line Psychological Experiments[END_REF]. Most displays or playbacks are written directly to a specific video or sound card to achieve the performance required of experimental software. Moreover, new video and sound cards are being released without consideration of DOS support, and more generally, recent PC computers no longer support DOS, which is quickly disappearing. It is clear that for a new application to have an acceptable life span, it has to be based on an up-to-date environment like Windows [START_REF] Forster | DMDX: A Windows Display Program with Millisecond Accuracy[END_REF].

Why are some solutions unable to obtain timing accuracy?

One way to limit timing inaccuracy on multitask operating systems is to give a high priority level to the dedicated program, which means that it is executed frequently in the computer task list. But even with this precaution, there is no guarantee that it will not be interrupted for a critical duration by another high-priority task (e.g. drivers managing hardware). The use of "timers" may be a solution for controlling timing accuracy. A timer is a way to execute a function regularly. For example, it can be used to detect whether the subject has responded (by systematic scanning of a response device). But standard timers running on Microsoft Windows have a resolution of about 50 ms, which is very poor for precisely measuring response time. Multimedia timers offer a better resolution (about 1 ms) but this precision is not guaranteed.

Why does command latency vary?

Another problem is related to the latency of some commands. The latency is the delay between the moment when a command is launched (e.g. changing a view, playing a sound) and the moment when the command is actually executed. For example, on a Pentium III 1.2 GHz, the C command PlaySound("stimulus.wav") to play back a wave file has a latency of 100 ms. This means that if a timer that measures response time is started just before the PlaySound instruction, it will have an error of about 100 ms because the actual audio stimulus started 100 ms later than planned. One solution could be to measure all the critical latencies once and for all, and take them into account in the results. But the problem is more complicated. In fact, these latencies can change, depending on a large set of factors, which means that this solution is inadequate.

How can timing accuracy be obtained?

First of all, the computer in charge of the perception test must be dedicated to this job only during the test [START_REF] Forster | DMDX: A Windows Display Program with Millisecond Accuracy[END_REF]. The operator must be sure to remove all tasks that can use up processing time (e.g. anti-virus, network, CD-Rom). But this is still not enough.

To obtain timing accuracy on the Perceval workstation, we use a well-suited technology called DirectX. DirectX provides a set of programming interfaces for designing high-speed applications while obtaining the benefits of direct access to the hardware. This software development kit includes several components. DirectSound provides low-latency mixing, hardware acceleration, and direct access to the sound device. DirectDraw allows one to directly manipulate display memory, hardware blitter, and the flipping surface support. These two components provide the capability of reducing the latency of playback and display functions at a high ratio because they are close to the hardware. Moreover, it is possible to prepare a set of time-consuming tasks before the real stimulation. For example, as a first step (preliminary), one can load a wave file or an image file in memory and then, as a second step (the real test), play back the sound and flip the screen image very rapidly. For comparison (see section 3.4), on a Pentium III 1.2 GHz, the DirectSound function "Play" has a latency of 0.5 ms. The third component used is DirectInput, which enables an application to gain access to input devices, including the mouse, keyboard, and joystick, even if the application is in the background. This capability can be implemented by using buffered data, which is a record of events that are stored until the application retrieves them.

MODULAR ARCHITECTURE

The Perceval package includes several specialized applications: stimulus creation tools, experiment design wizard, subject manager, test module, scoring module, etc.

Tests are controlled by scripts, a technique already used in the S.O.A.P. system [START_REF] Howard-Jones | SOAP' -A Speech Output Assessment Package for Controlled Multilingual Evaluation of Synthetic Speech[END_REF]. In fact, the system's core software resembles that of a lexical-syntactic analyzer, which reads and interprets the script files sent to it. The execution sequence (trial) is modified according to the commands and data received. This type of operation provides a great deal of flexibility and supports a wide variety of tests. A script is a text file composed of different sections: information, data, trial events, and settings (Figure 1). These files can be edited manually by an advanced user who knows the syntax. A novice can also follow the instructions given by an application design wizard with a menu-driven interface.

For each subject, a response file is generated containing the necessary information for processing the data. These files are written in table format to make it easy to export them to spreadsheets (Figure 2) or statistical analysis packages.

EXPERIMENTAL DESIGN AND RESULTS

Script example adapted to minimal pairs (Figure 1)

This test consists in playing the sounds, displaying word pairs, and recording the subject's responses. The results can be improved by including the features and contexts of the consonants tested. Usually, a script is divided into four major sections: The [TRIAL_EVENTS] section defines the event sequence for each trial played. The first order in this section is always BEGIN and last is always END. The command on line X20 displays the contents of the first column (#1) of the trial data on the screen. In the above example, if the trial played is TRIAL1, then 1)main 2)bain is displayed. The command on line X30 plays the wave file located in the second column (#2) of the trial data. In our example, if the trial played is TRIAL1, the sound bain.wav is played. The command on line X40 starts the recording of the response. In this example, two seconds are left for the subject to respond.

The [SETTINGS_GROUP1] section defines the experiment configuration (display, trial order, input, pause, etc.). Several groups can be included in the same script.

At the beginning of the experiment, instructions are displayed on the screen. The instructions are contained in a text file (*1).

A training phase can be executed. The order of the trials to be run is defined in *2. During the actual test phase, the order of the trials can be fixed or randomized (*3). The format (font, size, etc.) of the text displayed on screen is set in *4. The authorized keys for responding are listed in *5: a standard keyboard or a button box can be used for the test. Depending on the type of test, the concept of correct/incorrect response can be introduced (*6). In our example, this information is in the third column (#3) of the trial data. The pause value (*7) defines the time between trials. The results written in the response file are formatted in *8. Figure 2 shows an example of a response file obtained with the previous script (Figure 1). A: $SUBJECT specifies the subject code. B: $TRIAL defines the trial number. C: #1 is the contents of the first column of the trial data for that trial. In this example, it is the text displayed. D: #2 is the contents of the second column of the trial data for that trial. In this example, it is the sound played. E: #3 is the contents of the third column of the trial data for that trial. In this example, it is the correct response. F:$RESPONSE is the response given by the subject. G:$ERROR specifies whether the subject's response is correct or incorrect ("ok" if it correct, "err" if not).

H: #4 is the contents of the fourth column of the trial data for that trial. In this example, it is the consonant feature. I: #5 is the contents of the fifth column of the trial data for that trial. In this example, it is the vocalic context of the consonant. J: $RTIME is the subject's reaction time.

Example of a script for varying the volume of a sound (Figure 3)

The subject must determine whether the sound heard (aaa.wav) seems strong or weak (Fort/Faible). We vary the volume of the sound using the VOLUME order. In the example shown in Figure 3, if the trial executed is the first one, the sound is played as it was recorded. If the trial is the second one, the sound is played at 3 db less than it was recorded.

[TRIAL_DATA] TRIAL1=<0> TRIAL2=<-3> TRIAL3=<-6>

[TRIAL_EVENTS] X10=BEGIN X20=DISPLAY_TEXT <1)Fort 2)Faible> X30=PLAY_SOUND <aaa.wav><VOLUME #1> X40=GET_INPUT X50=END The subject listens to one of the two words ("bêle" or "bête") for a variable amount of time (gating). The beginning (TIME_BEGIN) and the end (TIME_END) of the sound listening time are defined. Here, the sound is played between the beginning and ending times defined by the third parameter in the trial data list.

[TRIAL_DATA] TRIAL1=<bêle><bele.wav><200> TRIAL2=<bêle><bele.wav><250> TRIAL3=<bêle><bele.wav><275> TRIAL4=<bête><bette.wav><200> TRIAL5=<bête><bette.wav><250> The subject sees images with an associated word, and must decide whether the word corresponds to the image. A sound is emitted to tell the subject whether the answer given was right or wrong (sound feedback).

[TRIAL_EVENTS] X10=BEGIN X20=DISPLAY_TEXT<1) bêle 2) bête> X30=PLAY_SOUND<#2><TIME_BEGIN 0><TIME_END #3> X40=GET_INPUT<DELAY 2000> X50=END
[TRIAL_DATA] TRIAL1=<catre><faute><catre.bmp> TRIAL2=<glace><faute><glace.bmp> TRIAL3=<horloche><faute><horloche.bmp> TRIAL4=<sourus><faute><sourus.bmp>

CONCLUSION

The computer-driven speech assessment system we developed is used to evaluate auditory-lexical decisions, phoneme monitoring, gating, phonetic categorization, word identification, voice quality, speech intelligibility, and so on. It was designed to automate the various operations involved in setting up experiments, presenting stimuli, and recording subject responses on a single PC computer. The experience acquired during the development process will be applied to enhancing the performance and reliability of our multiple-listener workstation, which can test as many as eight subjects at a time. PERCEVAL is available on the website of the "Speech and Language Laboratory" at the following address: www.lpl.univ-aix.fr

Figure 1 .

 1 Figure 1. Script exampleThe [INFORMATION] section contains general information: the name of the author, the creation date, the version used, the title given. The [TRIAL_DATA] section defines all data one needs during the experiment for each trial executed. Each line of data is divided into columns. The [TRIAL_EVENTS] section defines the event sequence for each trial played. The first order in this section is always BEGIN and last is always END. The command on line X20 displays the contents of the first column (#1) of the trial data on the screen. In the above example, if the trial played is TRIAL1, then 1)main 2)bain is displayed. The command on line X30 plays the wave file located in the second column (#2) of the trial data. In our example, if the trial played is TRIAL1, the sound bain.wav is played. The command on line X40 starts the recording of the response. In this example, two seconds are left for the subject to respond.

Figure 2 .

 2 Figure 2. Response file for minimal pairs (Excel file)

Figure 3 .

 3 Figure 3. Variation of a sound's volume 5.3 Example of a script using gating (Figure 4)

Figure 4 .

 4 Figure 4. Use of gating 5.4 Example of a script using feedback (Figure 5)

Figure 5 .

 5 Figure 5. Use of feedback


Alain Ghio
email: alain.ghio@lpl.univ-aix.fr

Bernard Teston
email: bernard.teston@lpl.univ-aix.fr

Evaluation of the acoustic and aerodynamic constraints of a pneumotachograph for speech and voice studies

To measure oral and nasal airflow during speech production, several conditions must be gathered. In this aim, we designed and built a pneumotachograph with particular care to optimise its response time, linearity and acoustical response. This flow meter is based on the grid flow meter principle with a small dead volume and specific linearization for the inhaled and exhaled airflow. A soft silicone rubber mask, pressed against the speaker's face prevents air leakage, without hindering articulatory movements. The acoustical distortions of the speech sound through the device are remedied by an adapted signal processing from its transfer function.

INTRODUCTION

The measurement of the aerodynamic parameters for the study of the articulatory mechanisms in the speech and voice production, particularly that of the oral and nasal flows examination, presents several problems which were never solved other than by compromises. To measure the air flow, a certain number of conditions must be joined together: 1-The flow meter must be bidirectional (inhaled and exhaled airflows). 2-Thermodynamic conditions of the air (temperature, composition, moisture) should not have effects to the measures. 3 -Turbulences with the lips (which depend on the type of phoneme) should not disturb measurements. 4-The dead volume of the sensor must be as small as possible to obtain a low acoustic distortion and a good response to fast variations of flow (especially for the oral air flow). 5-The flow meter should not have a resistance to the air flow likely to disturb the operation of the vocal tract. 6-The flow meter must be correctly adapted to the anatomy of the speaker to avoid leakages, source of errors of measurement, without obstructing the articulatory movements. 7-The device must be compatible with the hospital standards hygiene conditions.

It is very rare, in instrumentation, to meet a sensor which satisfies so many contradictory conditions. Some of them were already posed for measurements on breathing and quite naturally, phoneticians used sensors of respiratory flow (pneumotachographs) associated with anaesthesia masks, with not very reliable results, because these devices were very badly adapted to the dynamics of the speech production. In spite of this projection, very few systems were worked out for the articulatory studies. To equip the EVA2™ system (SQLab, www.lpl.univaix.fr/~sqlab), we produced an airflow sensor matching the main part of the conditions developed previously. The object of this study is to evaluate and correct the acoustic and aerodynamic distortions generated by this device.

PNEUMOTACHOGRAPHS

The pneumotachographs or PTG are sensors of air flow which measure the gaseous exchange in breathing. A good linearity was the main challenge for a long time. This quality was reached satisfactorily only by Fleisch [Fle25], whose PTG had a great dead volume and a significant resistance to the flow. To reduce this last disadvantage, Lilly [START_REF] Lilly | Flowmeter for recording respiratory flow of human subjects[END_REF] proposed the first PTG with grid whose various evolutions were used for many phonetic studies in particular by Guilbert [START_REF] Guilbert | Oral airflow during stop consonant production[END_REF]. Rothenberg [START_REF] Mr | Measurement of airflow in speech[END_REF] was the first of the rare attempts to carry out a PTG adapted to the particular constraints of aerodynamic measurements in speech production. With an aim of solving the problem of the response time of the PTG, he integrated the grid of measure to the walls of an anaesthesia mask. Always diffused, this system cannot really be regarded as a precise flow meter of measurement. Other principles of PTG were proposed such as those with hot wire or Doppler effect by ultrasounds [START_REF] Buess | Design and construction of has pulsed ultrasonic air flowmeter[END_REF], but they appeared incompatible with certain conditions developed in the introduction.

The aerophonometer

The PTG we study derives directly from the principle described by Teston [START_REF] Teston | Etude d'un aérophonomètre de grande dynamique et faible constante de temps[END_REF][START_REF] Teston | A new airflowmeter design for the investigation of speech production[END_REF]. Taking into account the specificity of its application and its realization, we baptized it "aerophonometer" to differentiate it from the standard PTG. This device is a PTG with grid (stainless steel wire of 200 µm diameter and a step of mesh of 250 µm), of size reduced (30 mms diameter and 20 mm length) to optimize its response time and its linearity in all the articulatory contexts. Thanks to its association with differential pressure transducers (Data Instrument DCXL) very sensitive and stable, it is able to apprehend a flow of the order of one cm3/s. The resistance of the grid is 10 Pa by dm3/s (litre/s), i.e. approximately 1% of the intra oral pressure of a normal subject, which does not disturb the operation of the vocal tract. Resistance was selected for a level of saturation of the sensor to the value of 10 dm3/s in forced breathing. These values represent a dynamics of 60 dB. To reduce the non-linearity of measurement caused by aerodynamic turbulences produced during speech production, the pressure tap is made in 8 points of the circumference of the measurement pipe and a grid of tranquillization (negligible in resistance) is laid out in front of the pressure taps. The sensor is made with synthetic material (Polyacetal) which has a very good mechanical resistance to sterilization and UV.

Adaptation to the anatomy of the speaker

This type of flow meter equips the system EVA2 ™ (Computerised Vocal Assessment, SQLab) which is a device designed for the recording and the study of many parameters of the speech and voice production: sound, pitch, intensity, airflows, pressures [START_REF] Teston | A diagnosis and rehabilitation aid workstation for speech and voice pathologies[END_REF]. The simultaneous recording of the sound and the airflows to the mouth and the nostrils requires the use of a "mouthpiece" which is a mechanical stand of the microphone, sensors of flow and flexible silicone mask allowing the sealing necessary to the reliability of measurements of oral air flow ( Fig 1). The nasal sensor is located under the oral sensor and vertically to follow the natural flow of the air to the nostrils. The flexible silicone masks, of a pleasant contact on the skin, take support on the speaker's face by releasing the nose and by taking the chin. They allow normal movements of the mandible in current speech, while preserving an excellent sealing. They are also easily sterilizable.

LINEARIZATION OF AIRFLOWS MEASUREMENTS

Position of the problem

The principle of measurement of the flow by a PTG consists in opposing a resistance to the air flow and to measure the pressure difference ∆P before and after. It In a general way, the relation between these values is formalized by: ∆P = R. D N .

The precise knowledge of the relation between ∆P and D makes it possible to bring back the measurement of flow to a measurement of differential pressure, operation technically easy to realize. Indeed, ∆P can be easily measured by differential pressure sensors. Those deliver a tension V proportional to ∆P. If the device is regulated in order to give a tension V calib for a flow gauged D calib , the measurement of the flow will be finally equal to:

(Eq.1) n calib calib n V D G V G D = = où
By applying a logarithmic transformation, we obtain:

(Eq.2) ( ) ( ) ( ) V Log n G Log D Log 1 + =

Measurements

The chain of instrumentation is composed of a blower which provides an adjustable and regular flow. This air flow is controlled by flow meters of high degree of accuracy TSI ™ model 4040 (for the flows from 0 to 5 dm 3 /s) or model 4140 (from 0 to 0.33 dm 3 /s). The air flow thus gauged is sent in the sensors whose output voltage is measured by a Fluke 45™ multimeter. A microcomputer PC controls the chain of instrumentation and records, via two ports series (RS232), the values provided by flow meter TSI ™ and the voltmeter. The oral flow sensor is not geometrically symmetrical because of the grid of tranquillization which has a significant effect. Also, measurements were carried out with blown and aspired air flows. This dissymmetry is marked for the sensor of nasal flow because of the air intake to the nostrils.

Results

A linear regression carried out to the logarithmic measures (cf Eq. 2) allows us to obtain the factors A and B of Eq.3. Thus, for the flow of expired oral air (range 2 dm 3 / S max), the measurement of the flow is obtained by the relation: D (in dm 3 / s) = 0.4436. V 0.654 where V is the tension in volts at exit of the sensors.

For values of low flow, the relation is linearised and tended towards a relation where D ( dm 3 / s) is locally proportional to the difference in pressure.

The adjustment of the function of correction for the flows inspired and expired with the mouth and the nostrils made it possible to obtain, in all the ranges of measurement, a linearity of the sensors corrected about 1%. The precision is validated by volumetric measures with a syringestandard of 3 liters (series 5530, Hans Rudolph Inc, the USA) for various rates of flow. Air volume is obtained by a temporal integration of the air flow signal.

Drift is another source of nonlinearity for such a sensitive sensor. The pressure sensors we use have a weak drift in temperature and it is useless to control it. However, an automatic reset is made before each recording to cancel this problem.

RESPONSE TIME OF THE FLOW METERS

Position of the problem

Response time of the flow meter is related to the phenomena which can delay the measure of a fast variation of flow. This delay is linked to the response time of the sensor or to its bandwidth in frequency. But the main cause of this delay is the dead volume of the flow sensor, great defect of the PTG based on respiratory measurement. The aerophonometer was built for that, with a geometrical volume of 8 cm3 and with very short pipes connections to the pressure sensors. These transducers (series DCXL, Data Instruments) have response times quite lower than 1 ms. The delay is also a function of the filtering applied to the signal of the flow before analog-to-digital conversion. Indeed, the fast and sensitive pressure sensors behave like poor microphones and can record the variations of acoustic pressures. To avoid aliasing, we use a 8 th order low pass Bessel filter with a cut-off frequency to 60 Hz. The indicial response of this function has no delay and bounce. So, it is well adapted to this use.

Measurements

To test the response of the flow meters, we produce an impulse of air flow thanks to a weight falling on an air bag. The noise of the shock is collected by the microphone of the mouthpiece and gives the precise moment of the beginning of the flow transient.

The total response time of the system is illustrated in figure 3 for the oral sensor. It is lower than the millisecond. The nasal sensor was optimized less in the

Position of the problem

The flow sensor, placed in the prolongation of the vocal tract, introduces a distortion of the speech sounds according to its own transfer function. Even if the measurement pipe is reduced to a minimum, it has a certain volume increased by the volume between the mask and the speaker's face. Additional cavities with the measurement pipe linked to the pressures sensors can also produce parasitic acoustic phenomena. Lastly, the device disturbs the natural acoustic radiation with the lips. In spite of the application of all the fundamental principles to avoid acoustic distortions as much as possible, there are always some which are necessary to correct.

Measurements

We measure the harmonic distortions of the mask-sensor thanks to an artificial head loudspeaker, designed for this use, under conditions identical to the use of the system ( It is noted that the mask-sensor unit behaves like a lowpass filter with a strong resonance towards 1.4 kHz caused by the additional cavity between the lips and the grid. The anti-resonance located at 3.5 kHz on the spectrum is probably the acoustic consequence of the side cavities of measurement pipe which collect each one 4 points of measurement of pressure. Their volume is about 1 cm3 and cannot be absorbed because of hygiene requirements (the sensor must be able to be sterilized easily). The corrections we do consist in applying, by digital filtering, an equalization whose frequential profile is the reverse of the transfer function of the system. The result of the correction is presented in Fig 6

CONCLUSION

This work of optimal devices adapted to speech production studies requires the realization of specific instruments of calibration as well as a constant rigour and a determination. Although satisfied of the general performances of the device, we consider in a future immediate to perfecting the whole system and, more particularly, the correction of the acoustic distortions.

Fig 1 .

 1 Fig 1. Speaker position with aerodynamic sensors.

  varies according to the speed of the fluid, i.e. its volumic flow D. For low values, the flow is laminar. The laws of the mechanics of the fluids show that, in this case, ∆P = R.D where R is resistance (equation of Poisseuille). For more significant values of flow, it becomes turbulent. The pressure loss ∆P is then proportional to the square of the flow (law of Venturi): ∆P = R. D².

  which translates a linear relation y = a.x+b (Eq.3) where y = Log(D), x=Log(V), a=1/n, B = Log(G) Our objective is to obtain empirically the parameters a and b from experimental measures ( Fig 2) .

Fig 2

 2 Fig 2 Pneumotachograph response curve of the expired oral air flow (for the range 2 dm 3 / s max)

  Fig 3a: Acoustic signal of the shock

  Fig 4). The signal source is provided by a white noise generator (Brüel & Kjaer model 1027). The recordings are carried out in soundproof room. Fig 5 represents the difference between the long-term spectrum of the response of the system and that of the artificial head.

Fig 4 Fig 5 Fig 6

 456 Fig 4 Artificial head loudspeaker to study acoustic distortion

  BIBLIOGRAPHY [Fle25] Fleisch A. (1925), "Der Pneumotachograph; ein Pageantry zur Geschwindigkeits-registrierung der Atemluft ", Pflugers Arch. Ges. Physiol., vol. 209, pp. 713-722


Fran De Aquino

Quantum Controller of Gravity

Keywords: Gravitation, Gravitational Mass, Inertial Mass, Gravity, Quantum Device

A new type of device for controlling gravity is here proposed. This is a quantum device because results from the behaviour of the matter and energy at subatomic length scale (10 m).

-20

From the technical point of view this device is easy to build, and can be used to develop several devices for controlling gravity.

Introduction

Some years ago I wrote a paper [START_REF] De Aquino | Mathematical Foundations of the Relativistic Theory of Quantum Gravity[END_REF] where a correlation between gravitational mass and inertial mass was obtained. In the paper I pointed out that the relationship between gravitational mass, , and rest inertial mass, , is given by g m 0 i m ( ) is the variation in the particle's kinetic momentum;

1 1 1 2 1 1 1 2 1 1 1 2 1 2 2 2 2 0 2 0 0 ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + - = = ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + - = = ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ Δ + - = =
is the electromagnetic energy absorbed or emitted by the particle;

is the index of refraction of the particle; W is the density of energy on the particle ; p Δ or that g g χ = ′

In the last years, based on these concepts, I have proposed some types of devices for controlling gravity. Here, I describe a device, which acts controlling the electric field in the Matter at subatomic level ( )

m x 20 10 - ≅ Δ
. This Quantum Controller of Gravity is easy to build and can be used in order to test the correlation between gravitational mass and inertial mass previously mentioned.

The Device

Consider a spherical capacitor, as shown in Fig. 1. The external radius of the inner spherical shell is , and the internal radius of the outer spherical shell is . Between the inner shell and the outer shell there is a dielectric with electric permittivity . The inner shell works as an inductor, in such way that, when it is charged with an electric charge q + , and the outer shell is connected to the ground, then the outer shell acquires a electric charge , which is uniformly distributed at the external surface of the outer shell, while the electric charge is uniformly distributed at the external surface of the inner shell (See Halliday, D. and Resnick, R., Physics, Vol. II, Chapter 28 (Gauss law), Paragraph 28.4). q q +

+ + + + V +q -q V 2 =0 V 1 =V r r b r a R 10
Under these conditions, the electric field between the shells is given by the vectorial sum of the electric fields Note that the electrostatic force, , between and will move the negative electric charges in the direction of the positive electric charges. This causes a displacement, , of the electric field,

F r q - q + x Δ b E r
, into the outer shell (See Fig. 2 (b)). Thus, in the region with thickness x Δ the intensity of the electric field is not null but equal to .

b E

The negative electric charges are accelerated with an acceleration, a , in the direction of the positive charges, in such way that they acquire a velocity, given by r

x a v Δ = 2 (drift velocity).
The drift velocity is given by [START_REF] Griffiths | Introduction to Electrodynamics[END_REF] ( )

2 2 2 nSe X R V nSe Z V nSe i v C + = = =
where is the positive potential applied on the inner shell (See Fig. 1 , into the outer shell. Thus, in the region with thickness x Δ the intensity of the electric field 2) can be rewritten in the following form:

); V fC X C π 2 1 =
is not null but equal to b E . r a + - r b E b E a E a E a E b E b E R =E a + E b E R = 0 E R = 0 E= 0 E= 0 ε V r a + - r b E b E a E a E a E b E b E R =E a + E b E R = 0 E R = 0 E= 0 E= 0 ε V Δx F If
( )

3 10 8 . 6 20 10 V nSe R V nSe i v - × = ≅ =
The maximum size of an electron has been estimated by several authors [START_REF] Dehmelt | A Single Atomic Particle Forever Floating at Rest in Free Space: New Value for Electron Radius[END_REF]4,[START_REF] Macken | Spacetime Based Foundation of Quantum Mechanics and General Relativity[END_REF]. The conclusion is that the electron must have a physical radius smaller than 10 -22 m * .

Assuming that, under the action of the force (produced by a pulsed voltage waveform,V ), the electrons would fluctuate about their initial positions with the amplitude

of F r m x 20 10 - 1× ≅ Δ
(See Fig. 3), then we get ( )

4 294 . 0 2 2 V v x a x t ≅ Δ = Δ = Δ However, we have that t T f Δ = Δ = 2 1 1
. Thus, we get ( )

5 7 . 1 V f = Now consider Eq. (1)
. The instantaneous values of the density of electromagnetic energy in an electromagnetic field can be deduced from Maxwell's equations and has the following expression ( )

m x 20 10 - ≅ Δ . t Δ E b t Δ E b m x 20 10 1 - × ≅ Δ E b E b E b - - - + + + F F - Equation (7) shows that v r = κ ω . Thus, v k B E r = = ω , i.e., H v vB E μ = =
Then, Eq. ( 6) can be rewritten in the following form:

( ) ( ) 9 2 2 1 2 2 2 1 H H v W μ μ μ ε + = For ωε σ << , Eq. (7) reduces to r r c v μ ε = Then, Eq. (9) gives 2 2 2 1 2 2 2 1 H H H c W r r μ μ μ μ μ ε ε = + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ =
This equation can be rewritten in the following forms:

( )

10 2 μ B W = or ( ) 11 2 E W ε = For ωε σ >>
, Eq. ( 7) gives ( )

12 2 μσ ω = v
Then, from Eq. (9) we get ( )

13 2 2 2 1 2 2 1 2 2 2 1 2 2 1 H H H H H W μ μ μ σ ωε μ μ μ μσ ω ε ≅ ≅ + ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ = + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = Since H v vB E μ = =
, we can rewrite (13) in the following forms:

( ) 14 2 2 μ B W ≅ or ( ) 15 4 2 E W ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ≅ ω σ
Substitution of Eq. (15) into Eq. ( 2), gives ( )

16 1 10 758 . 1 1 2 1 1 256 1 2 1 1 4 4 1 2 1 0 4 3 2 3 27 0 4 3 2 3 2 3 0 0 2 4 3 2 i r i r i g m E f m E f c m E f c m ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ × + - = = ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + - = = ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + - = - ρ σ μ ρ σ μ π μ ρ π σ μ
Using this equation we can then calculate the gravitational mass, , of the region with thickness ( ) ( )

x g m Δ x Δ ,
19 1 10 3 . 1 1 2 1 0 3 4 2 x i x g m f V m Δ - Δ ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - × + - =
Equation [START_REF] Macken | Spacetime Based Foundation of Quantum Mechanics and General Relativity[END_REF] shows that there is a correlation between V and to be obeyed, i.e.,

f V f 7 . 1 =
. By substituting this expression into Eq. ( 19), we get , the gravity will be .

x Δ g 2 χ † Note that the frequency must be greater than 1Hz

(See text above Eq. ( 3)). (Eq. ( 5)), then it is necessary to put a synchronizer before the pulse generator (See Fig. 5), in order to synchronize V with f .Thus, when we increase the voltage, the frequency is simultaneously increased at the same proportion, according to Eq. ( 5).

Pulse Generator V, f

Synchronizer V 0-1. 2kV 0-2. 04kHz g Fig. 4 -Experimental Set-up using a Quantum Controller of Gravity (QCG).

  shown that, if the weight of a particle in a side of a lamina is lamina) then the weight of the same particle, in the other side of the lamina is respectively, the gravitational mass and the inertial mass of the lamina). Only when equal in both sides of the lamina. The lamina works as a Gravity Controller.

Fig. 1 -

 1 Fig.1 -Spherical Capacitor -A Device for Controlling Gravity developed starting from a Spherical Capacitor.

.

  by the inner shell and the outer shell. Since they have the same direction in this region, then one can easily show that the resultant intensity of the electric field for In the nucleus of the capacitor and out of it, the resultant electric field is null because a E r and b E r have opposite directions (See Fig. 2(a)).

Fig. 2 -

 2 Fig.2 -The displacement, x Δ , of the electric field, b E r

  the shells are made with Aluminum, with the following characteristics: ,

;Fig. 3 -

 3 Fig.3 -Controlling the Electric Field in the Matter at subatomic level ( )

  in the outer shell. We have already seen that the electric field in this region is b capacitor; V is the potential applied on the inner shell (See Fig.1 and 3). Thus, Eq. (17) can be rewritten as follows

Fig. 4 -χ

 4 Fig.4 -The shell with thickness x Δ works as a Quantum Controller of Gravity.


Channa Mey

Artem Metlin

Veasna Duong

Sivuth Ong

In Sotheary

Paul F Horwood

Jean-Marc Reynes

Hervé Bourhy

Arnaud Tarantola

Philippe Buchy
email: buchyphilippe@hotmail.com

Evidence of two distinct phylogenetic lineages of dog rabies virus circulating in Cambodia SHORT TITLE: Phylogeny of rabies virus in Cambodia CONTENTS CATEGORY: Standard (full-length)

Keywords: rabies virus evolution, rabies epidemiology, canine rabies, Cambodia, South-East Asia

This first extensive retrospective study of the molecular epidemiology of dog rabies in Cambodia included 149 rabies virus (RABV) entire nucleoprotein sequences obtained from 1998-2011. The sequences were analyzed in conjunction with RABVs from other Asian countries. Phylogenetic reconstruction confirmed the South-East Asian phylogenetic clade comprising viruses from Cambodia, Vietnam, Thailand, Laos and Myanmar. The present study represents the first attempt to classify the phylogenetic lineages inside this clade, resulting in the confirmation that all the Cambodian viruses belonged to the South-East Asian (SEA) clade. Three distinct phylogenetic lineages in the region were established with the majority of viruses from Cambodia closely related to viruses from Thailand, Laos and Vietnam, forming the geographically widespread phylogenetic lineage SEA1. A South-East Asian lineage SEA2 comprised two viruses from Cambodia was identified, which shared a common ancestor with RABVs originating from Laos. Viruses from Myanmar formed separate phylogenetic lineages within the major SEA clade. Bayesian molecular clock analysis suggested that the time to most recent common ancestor (TMRCA) of all Cambodian RABVs dated to around 1950. The TMRCA of the Cambodian SEA1 lineage was around 1964 and that of the SEA2 lineage was around 1953. The results identified three phylogenetically distinct and geographically separated lineages inside the earlier identified major SEA clade, covering at least five countries in the region. A greater understanding of the molecular epidemiology of rabies in South-East Asia is an important step to monitor progress on the efforts to control canine rabies in the region.

Highlights:

We performed the first extensive rabies virus evolutionary study in Cambodia;

We generated and analyzed 149 Cambodian rabies virus nucleoprotein sequences;

We identified 3 phylogenetic lineages in the major South-East Asian clade.

Introduction

The developing world is severely affected by rabies, causing deaths in humans as well as in animals and resulting in significant economic losses. Indeed, it is estimated that globally canine rabies causes approximately 59,000 human deaths, over 3.7 million disability-adjusted life years (DALYs) and 8.6 billion USD in economic losses annually [START_REF] Hampson | Global Alliance for Rabies Control Partners for Rabies Prevention[END_REF].

Rabies virus (RABV) belongs to the family Rhabdoviridae of the order Mononegavirales which comprises at least 10 genera (http://www.ictvonline.org). The Lyssavirus genus includes classical RABVs, found in a wide range of different animal species throughout the world. All 13 other species of Lyssavirus except two have been isolated from bats, originating from different locations in Europe, Asia, Africa and Australia. RABV has a negative-sense RNA genome of about 12 kb encoding five proteins [START_REF] Tordo | Walking along the rabies genome: is the large G-L intergenic region a remnant gene?[END_REF][START_REF] Tordo | The rabies virus genome: an overview[END_REF]. The nucleoprotein gene (N-gene) is widely used for molecular characterization and phylogenetic analysis of RABVs [START_REF] Kissi | Genetic polymorphism in the rabies virus nucleoprotein gene[END_REF]. The large number of N-gene sequences published makes this an obvious region of the genome to study [START_REF] Johnson | Phylogenetic comparison of the genus Lyssavirus using distal coding sequences of glycoprotein and nucleoprotein genes[END_REF][START_REF] Chupin | Genetic characterization of the rabies virus field isolates detected in Russian Federation within the period 2008-2011[END_REF].

The N-gene is also suitable for estimations of evolutionary rate and subsequent analysis of evolutionary history [START_REF] Holmes | Genetic constraints and the adaptive evolution of rabies virus in nature[END_REF][START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF].

Most South-East Asian countries are affected by rabies, but with the exception of a few countries, there is limited accurate epidemiological information available. According to the World Organization for Animal Health (OIE), there were 15 animal rabies outbreaks recorded in Vietnam in 2012 with 86 rabies cases found in dogs in five different administrative divisions of the country. For Laos it was only mentioned that rabies was present in 2012 with no indication of the number of outbreaks/cases. In Myanmar, the number of human rabies cases per year was estimated at about 1000 while in Thailand it was less than 25 [START_REF] Gongal | Human rabies in the WHO Southeast Asian Region: Forward Steps for Elimination Advances in Preventive Medicine[END_REF].

The only reliable data on human and animal rabies in Cambodia are available from the Institut Pasteur in Cambodia (IPC) where rabies diagnosis activity has been performed routinely since 1998, but laboratory-confirmed cases of rabies in dogs were recorded since the 1970s [START_REF] Reynes | Characterization and observation of animals responsible for rabies post-exposure treatment in Phnom Penh, Cambodia[END_REF]. In the second half of the 1970s, dogs which are the most important vector for human rabies in Asia [START_REF] Edited Bymeslin | Executive summary[END_REF] almost disappeared from Cambodia because of starvation. However, dog numbers rebounded in subsequent decades and rabies consequently became a serious public health concern. During the mid-1980s, the number of patients who have sought medical assistance after a dog bite has been at least 4000 per year, but the number of reported rabies deaths has been very low, partly because of improved postexposure treatment delivery. Deaths occurring at home are usually not reported to health services. From 1982 to 1991 a total of 5437 animal bites were recorded and 51 human rabies deaths were reported (WHO/Rab. Res./93.44). Since 1994, rabies stopped being a notifiable disease in Cambodia. Thus, the only data available after 1994 is from the IPC, which established a rabies post-exposure treatment centre in 1995 [START_REF] Reynes | Characterization and observation of animals responsible for rabies post-exposure treatment in Phnom Penh, Cambodia[END_REF]. From mid-1995 to 2007, 149,224 post-exposure prophylaxis (PEP) treatments were provided at the IPC Rabies Prevention Centre [START_REF] Reynes | Characterization and observation of animals responsible for rabies post-exposure treatment in Phnom Penh, Cambodia[END_REF]Ly, 2009; WHO/EMC/ZOO/96.8) and still 67 fatal human cases were reported countrywide following dog bites.

During 1998-2007, 610 animal samples (49% of all specimens tested) originating mostly from dogs with suspected rabies were found rabies-positive [START_REF] Ly | Rabies situation in Cambodia[END_REF]. A predictive model established by [START_REF] Ly | Rabies situation in Cambodia[END_REF] estimated that 810 human rabies cases occurred in 2007 alone. In 2008-2013, rabies was recorded mostly in dogs (n=932, 99.3%). Rabies cases in cats (n=3, 0.3%), bovines (n=3, 0.3%) and swine (n=1, 0.1%) were detected in very rare occasions (Institut Pasteur in Cambodia, unpublished data).

At present, the dog population in Cambodia is estimated to exceed 5 million animals and nearly one-third of the 250 Cambodian households surveyed in a preliminary study could recall at least one dog bite experience during 2004-2009 (Institut Pasteur in Cambodia, unpublished). Overall, the risk proportion (total number dog bites/total human population) in Cambodia was 5.6% (75/1339). This equates to an annual incidence of 1120 dog bites per 100 000 people [START_REF] Lunney | Knowledge, attitudes and practices of rabies prevention and dog bite injuries in urban and peri-urban provinces in Cambodia[END_REF]. Other mammals like mongoose, bats or other wild animals are sometimes suspected in Asia to contribute to maintaining the circulation of rabies. In Cambodia, [START_REF] Reynes | Serologic evidence of Lyssavirus infection in bats, Cambodia[END_REF] found serological evidence of Lyssavirus infection in bats. Several studies on RABV molecular epidemiology in Asia [START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF], namely in Vietnam [START_REF] Nguyen | Molecular epidemiology of rabies virus in Vietnam (2006-2009)[END_REF][START_REF] Yamagata | Molecular epidemiology of rabies in Vietnam[END_REF], Thailand [START_REF] Ito | Molecular epidemiology of rabies in Thailand[END_REF] and Laos [START_REF] Ahmed | Molecular epidemiology of rabies viruses circulating in two rabies endemic provinces of Laos, 2011-2012: regional diversity in Southeast Asia[END_REF] were conducted previously.

The objective of the present study was to perform a retrospective molecular and phylogenetic characterization of RABVs from dogs originating from different parts of Cambodia, to compare them to strains originating from other South-East Asian countries and to analyze the circulation of dog RABVs in this region.

Materials and Methods

2.1. Samples and Viruses. Dog brain samples received by the IPC Virology Unit from different regions of Cambodia between 1998 and 2011 were included in the study. In total, 149 rabies-positive samples obtained from dogs originating from 20 administrative regions of Cambodia were selected by location and year of sampling (Supplementary Table 1). Dog's heads were usually referred to the IPC lab from people who were consulting for PEP following the animal bite. The animal samples were tested by a standard direct fluorescent antibody test (FAT) [START_REF] Dean | The fluorescent antibody test[END_REF] using an anti-rabies nucleocapsid conjugate (Anti-Rabies Nucleocapsid Conjugate; lyophilizied, adsorbed #357-2112, Bio-Rad, Marnesla-Coquette, France) according to the manufacturer's instructions.

2.2. RNA purification, Reverse Transcriptase PCR and Nucleotide Sequencing. The diagnostic method for the initial detection of the RABV from clinical samples used a hemi-nested RT-PCR targeting a conserved region of the L-gene to confirm results of FAT as previously published [START_REF] Dacheux | A reliable diagnosis of human rabies based on analysis of skin biopsy specimens[END_REF]. The target size of amplicons obtained after the second round of PCR was approximately 250 base pairs. Negative and positive control samples were included successively for each extraction, reverse transcription, and for all PCR steps.

To generate amplicons for sequence analysis, two overlapping RT-PCRs were conducted targeting the N-gene [START_REF] Yang | Molecular characterization of Korean rabies virus isolates[END_REF], which has been extensively used for phylogenetic analyses of RABVs. The two sequencing PCRs (Table 1) were conducted using the same reaction conditions and cycling protocols as previously described [START_REF] Dacheux | A reliable diagnosis of human rabies based on analysis of skin biopsy specimens[END_REF]. The first PCR, using the primers RVN-

71F

(5'-ATGGATGCCGACAAGATTGTATTC- Phylogenetic analyses were conducted on RABV sequences including selected sequences sourced from the GenBank Database. Codon positions included were 1st+2nd+3rd+Noncoding. All positions containing gaps and missing data were eliminated. There were a total of 1346 positions in the final dataset for entire N-gene sequences, and 381 positions in the final dataset for partial N-gene sequences.

Evolutionary analyses were conducted in MEGA6 [START_REF] Tamura | MEGA6: Molecular Evolutionary Genetics Analysis version 6.0[END_REF].

The evolutionary history of the entire N-gene sequences was inferred using the Maximum Likelihood method. The evolutionary distances were computed using the Tamura 3-parameter model [START_REF] Tamura | Estimation of the number of nucleotide substitutions when there are strong transition-transversion and G + C-content biases[END_REF]. The percentage of trees in which the associated taxa clustered together is shown next to the branches. Initial tree(s) for the heuristic search were obtained by applying the Neighbor-Joining method to a matrix of pairwise distances estimated using the Maximum Composite Likelihood (MCL) approach. A discrete Gamma distribution was used to model evolutionary rate differences among sites (5 categories (+G, parameter = 0.3114)). The tree is drawn to scale, with branch lengths measured in the number of substitutions per site. The bootstrap test of phylogeny with 1000 replicates was performed.

The analysis involved 192 nucleotide sequences, obtained in the present study (Supplementary Table 1)

and retrieved from public databases (Supplementary Table 2).

For the partial N-gene, sequences were analyzed as above with some modifications. The evolutionary history was inferred using the Minimum Evolution (ME) method [START_REF] Rzhetsky | A simple method for estimating and testing minimum evolution trees[END_REF].

The tree was drawn to scale, with branch lengths in the same units as those of the evolutionary distances used to infer the phylogenetic tree. The evolutionary distances were computed using the number of differences method [START_REF] Nei | Molecular Evolution and Phylogenetics[END_REF] and are in the units of the number of base differences per sequence. The ME tree was searched using the Close-Neighbor-Interchange (CNI) algorithm [START_REF] Nei | Molecular Evolution and Phylogenetics[END_REF] at a search level of 1. The Neighbor-Joining algorithm [START_REF] Saitou | The neighbor-joining method: A new method for reconstructing phylogenetic trees[END_REF] was used to generate the initial tree. Bayesian Skyline Plot (BSP) was used as a coalescent prior. This provides a piecewise graphical depiction of changes of relative genetic diversity through time (N e τ), where N e is the effective population size and τ is the host-to-host generation time. We also utilized both strict and relaxed (uncorrelated lognormal) molecular clocks. A previous analysis of RABV evolution found that the uncorrelated lognormal relaxed clock and constant population size model provided a better fit to the data [START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF], and this was confirmed in a preliminary analysis of our data, although the parameter values estimated were similar under both models. The degree of uncertainty in each parameter estimate is provided by the 95% highest posterior density (HPD) values, while posterior probability values provide an assessment of the degree of support for each node on the tree. All Markov Chain Monte Carlo chains were run for sufficient time at 200 million steps to ensure statistical convergence, with 10% removed as burn-in.

Results

Sampling.

The geographical origins of the Cambodian dog samples analyzed in this study are shown in Figure 1. As there was no information on the exact location where the bite occurred, the cases were attached to the local district health center, which in Cambodia is usually located no more than 5 kilometers from the furthest point of the district.

Phylogenetic analyses.

The complete N-gene sequences of 149 RABV isolates collected in 20 provinces of Cambodia were generated and submitted to GenBank (Supplementary Table 1).

The Model Selection analysis with MEGA6 package revealed that for the dataset of entire and partial nucleoprotein sequences included in the present study, the Tamura 3 -parameter model with a discrete Gamma distribution (+G) had the lowest Bayesian Information Criterion (22395.836 for entire and 10845.81921 for partial nucleoprotein sequences) and was considered to describe the substitution pattern the best.

The complete N-gene sequences (1492 nucleotides in length, positions 71 to 1562, cut with the GeneDoc package to the length of 1346 nt to fit the equal size of the sequences of the reference strains, taken from the public databases) were used to build a Maximum Likelihood tree (Fig. 2). Phylogenetic analysis was then conducted on 194 partial N-gene sequences generated in the present study and other strains originating from South-East Asian countries (Supplementary Table 2).

Phylogenetic reconstruction identified the

Rabies virus 9337SWI isolated from a bat in 1993 (EBLV2) was used as an out-group to build this tree.

This analysis allowed a comparison with the shorter sequences obtained by Nguyen et al.

(2011) from

Vietnamese RABVs, resulting in further support for the clade and lineage groupings (Supplementary Fig.

Evolutionary analysis.

Bayesian molecular clock analysis performed on complete N-gene sequences suggested that the TMRCA of all 155 Cambodian rabies viruses dated back to 1950 (95% HPDs: 1935 to 1961; Fig. 3). The TMRCA of the Cambodian SEA1 lineage was around 1964 (95% HPDs: 1955 to 1971) and that of SEA2 lineage was around 1953 (95% HPDs: 1938 to 1969). The nucleotide substitution rate was estimated at 3.13 ×10 -4 substitutions/site/year (s/s/y; 95% HPD: 2.46 × 10 -4 -3.79 × 10 -4 ). Using Bayesian skyline plots as an inference tool, the patterns of relative genetic diversity in Cambodian rabies virus increased to its maximum approximately in 2005 followed by a small decline afterward.

Discussion

Domestic/stray dog rabies is still widespread throughout the world, occurring in over 80 countries/territories (predominantly in the developing world) and is implicated in more than 99% of all cases of human rabies (WHO, 2012). The present study aimed to provide greater clarity on dog rabies molecular epidemiology in South-East Asian countries using the sequences generated here in addition to the limited number of sequence data publicly available. More specific attention was paid to Cambodia, which is severely affected by rabies but where the epidemiology of the disease remains unclear.

Currently, there is no active rabies surveillance system in Cambodia, and the existing data is only available through the rabies diagnosis efforts and the PEP treatment activities of the IPC. As such, reliable analysis of the precise dynamics of dog rabies in Cambodia by year cannot be performed, but the available data confirms that the disease exists and is widespread throughout the country. summarized that RABVs from Thailand, Malaysia, Cambodia, Laos and Vietnam originated from a common ancestor; and the lineages circulating in Laos were closely related to each other but were not the same as in other South-East Asian countries.

In this first extensive study of RABV molecular epidemiology in Cambodia, a large pool of dog RABV nucleoprotein gene sequences were generated from strains detected in 20 provinces of Cambodia between 1998 and 2011. Figure 1 shows that these viruses were nevertheless mostly collected in the capital and its surrounding provinces. This is because the only rabies diagnostic facility -the Rabies National Reference Centre (IPC, Virology Unit) -and the major post-exposure treatment facility (IPC, Rabies Vaccination Centre) are both located in Phnom Penh. As a result, the samples collected from suspected rabid dogs were rarely received from the districts that are very distant from the capital.

Nevertheless, IPC's Rabies Reference Centre adequately covers the South-East region of Cambodia, which includes approximately 60-70% of Cambodian's human and canine populations. Cambodian viruses were included for phylogenetic analysis with RABV sequences from China, Laos, Myanmar, Thailand and Vietnam. Viruses originating from dogs were given preference to be included in analyses.

Some human RABV sequences from Vietnam were also included because of the lack of suitable representative sequences from dog origin from this country. However, as dogs are the major rabies reservoir in Asia [START_REF] Edited Bymeslin | Executive summary[END_REF], and most human cases had a history of a dog bite [START_REF] Nguyen | Molecular epidemiology of rabies virus in Vietnam (2006-2009)[END_REF], this should not have introduced a significant bias into the analysis.

The results presented in this study show that the major SEA phylogenetic clade covers at least of rabies in South-East Asia, further studies need to be conducted on an international level. To better define the molecular characteristics of the SEA clade and achieve a more unified terminology, more

Asian RABV sequences should be generated, eventually through a regional initiative. These analyses should also include rabies viruses from wildlife and bat Lyssaviruses of South-East Asian origin, the existence of which was suggested based on serological evidence [START_REF] Reynes | Serologic evidence of Lyssavirus infection in bats, Cambodia[END_REF].

The wide distribution of the SEA1 lineage across SE Asia suggests a relatively frequent mixing of dog populations throughout the region. This may provide new insights to assist the human and animal health authorities in these countries to take suitable actions to reduce dog circulation in their efforts to eradicate rabies. The Asian Rabies Elimination Strategy (ASEAN, 2015) has outlined the ambitious target of eliminating canine rabies in South-East Asia by 2020. Further clarity on the circulation of RABVs and dog populations will be important information to inform dog population management strategies.

2. 3 .

 3 Evolutionary analysis. To determine the comparative evolutionary dynamics of RABVs from Cambodia, we used the BEAST package (Drummond and Rambaut, 2007) to estimate the rate of nucleotide substitution per site, the Time to Most Recent Common Ancestor (TMRCA), and the Maximum Clade Credibility (MCC) tree of the complete coding region of N-gene sequences (Drummond and Rambaut, 2007) utilizing information on the exact year of sample collection. This analysis used the GTR model of nucleotide substitution as well as a different substitution rate for each codon position. A

  South-East Asian (SEA) clade that included all the nucleoprotein gene sequences of RABVs from Cambodia, Vietnam, Thailand, Laos and Myanmar. Inside the SEA clade 3 distinct phylogenetic lineages (SEA1, SEA2 and SEA3) were distinguished based on the bootstrap values obtained for the ancestral nodes of each cluster. Most of the viruses from Cambodia were closely related to viruses from Thailand, Vietnam and Laos and formed the lineage SEA1. We identified three putative sub-lineages within the SEA1 lineage, and two of them were geographically restricted to Thailand and Laos (putative sub-lineages 2 and 3, respectively). The SEA2 lineage was comprised of two viruses from Cambodia and a group of RABVs from Laos. One of these Cambodian RABVs was isolated from a dog in Phnom Penh in 2005; the second was isolated from a dog in 1997 from Banlung, the capital city of the Ratanakiri province in the North-East part of Cambodia, a region that shares common borders with Laos and Vietnam. This strain was already previously described by[START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF]. The two Cambodian SEA2 strains shared a common ancestor with RABVs isolated in Laos in 1999-2002[START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF] and in 2011-2012[START_REF] Ahmed | Molecular epidemiology of rabies viruses circulating in two rabies endemic provinces of Laos, 2011-2012: regional diversity in Southeast Asia[END_REF]. Viruses from Myanmar formed a separate phylogenetic lineage SEA3.

1 )

 1 . The SEA1 lineage, which was initially proposed by Yamagata et al. (2007), was comprised of viruses from Cambodia, Thailand, Laos and Vietnam and corresponded to the subgroup 1a suggested by Nguyen et al. (2011); whereas the SEA2 and SEA3 lineages established in the present study were distinct from the subgroup 1a and 1b, but still contained within the broader Group 1 described previously by Nguyen et al. (2011).

  Previously a large-scale research study covering South-East Asian countries provided the first definitive evidence for a widely distributed Asian clade comprising viruses from Myanmar, Thailand, Laos, Cambodia, Vietnam and some other countries, but the lineages in which the strains were clustering could not be determined[START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF].Ito et al. (1999) performed the first detailed study of rabies molecular epidemiology in Thailand in which six lineages of RABV were identified. At least four of them were endemic to Central (T1C), Southern (T2S) and Northern (T1N1 and T1N2) areas of Thailand whereas the two other lineages (T1CS and T2NC) overlapped different areas. Yamagata et al. (2007) proposed the SEA1 lineage of RABV including viruses from Thailand and Vietnam. Later Nguyen et al. (2011) established the existence of two major RABV lineages in Vietnam (Groups 1 and 2), with Group 1 further subdivided into two subgroups (1A and 1B). Subgroup 1A represented the SEA1 lineage previously described by Yamagata et al. (2007). Subgroup 1B was comprised of viruses from Vietnam, China and Philippines; and Group 2 included viruses from China and Vietnam. Ahmed et al. (2015)

  five countries in the South-East Asian region, including Cambodia, Vietnam, Thailand, Laos and Myanmar. This data, together with the timescale of RABV evolution, are in agreement with the results published earlier by[START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF]. In the present study, we attempted to investigate the presence of distinct phylogenetic lineages inside this large clade.Most of the lineages identified inside the SEA clade were associated with their geographical origins. All of the viruses from Cambodia, except two (EU086170:[START_REF] Bourhy | The origin and phylogeography of dog rabies virus[END_REF] and P0709509: this study), clustered with RABVs from Thailand, Vietnam and Laos within the SEA1 lineage, which appears to be the most widespread lineage throughout South-East Asia. We propose the existence of at least 3 putative sub-lineages inside the SEA1 lineage. The sub-lineage 1 comprised RABVs from Cambodia, Thailand, Vietnam and Laos. Two dog RABVs (GQ303556 & JN786877) clearly formed a distinct sub-lineage 2. Failure of rabies PEP associated with unusual manifestation of the disease was previously described with the strain GQ303556[START_REF] Shantavasinkul | Failure of rabies postexposure prophylaxis in patients presenting with unusual manifestations[END_REF]. The second virus, JN786877, after a series of passages in BHK cells rapidly killed mice inoculated by intramuscular route[START_REF] Virojanapirom | Molecular analysis of the mutational effects of Thai street rabies virus with increased virulence in mice after passages in the BHK cell line[END_REF]. Although only a small number of viruses from this phylogroup have been described, the potential association of these strains with an increased virulence should be further explored. Sub-lineage 3 comprised RABVs isolated in Laos (2011-2012) described by[START_REF] Ahmed | Molecular epidemiology of rabies viruses circulating in two rabies endemic provinces of Laos, 2011-2012: regional diversity in Southeast Asia[END_REF]. The SEA1 lineage correlates with the group 1/subgroup 1a, determined earlier by Nguyen et al.(2011) and the South-East Asian lineage 1 identified by[START_REF] Yamagata | Molecular epidemiology of rabies in Vietnam[END_REF].The SEA2 lineage seemed to be distributed in Cambodia and Laos, whereas the distribution of the SEA3 lineage appeared to be limited to Myanmar. The newly established SEA2 lineage comprised several RABVs from Laos and only two Cambodian strains: one from Phnom Penh and one from Banlung, located in the North-East of Cambodia near the Laos border. Taking into account that Laos RABVs included in the SEA2 lineage were either from the Vientiane capital, or from the Champasak province which is bordering Cambodia, the circulation of this lineage in Northern Cambodia and Southern Laos with occasional spread to the capitals of the two countries seems plausible. Further studies need to be conducted to confirm the precise origin and distribution of this lineage as there is a clear bias of recruitment of samples towards the capital cities. Only few strains originating from provinces bordering Laos and Vietnam were available for our study: one from the province of Preah Vihear, detected in 2009, and one from Banlung, Ratanakiri province, detected in 2010. It was unfortunately not possible to find any strains from the northern provinces of Cambodia, close to Laos, and thus interpretations must be exercised with caution. In Myanmar, only one lineage of RABVs was identified, comprising the SEA3 group. However, only three RABV sequences from Myanmar were available for analysis, all originating from 1999, so additional studies with more RABV sequences from this country are needed to confirm this observation. The SEA1 lineage was previously identified by Yamagata et al. (2007) although the bootstrap values were very low (32-40). It comprised viruses from Thailand and Vietnam and the authors proposed that because Cambodia and Laos are situated between Vietnam and Thailand, there was a possibility that RABVs from Cambodia and Laos were genetically related to the SEA1 lineage. This hypothesis is now confirmed by the present study. Such a wide distribution suggests the presence of suitable conditions allowing the intensive circulation of the viruses belonging to that lineage in the dog populations of Thailand, Laos, Cambodia and Vietnam. Illegal trade of stray dogs, essentially for meat consumption, between these four countries could well have contributed to the regional spread of RABVs from the SEA1 lineage (http://www.animals24-7.org/2014/06/08/dog-meat-traffic-still-spreads-rabies-invietnam/; Ahmed et al., 2015). The evolutionary analysis of the data generated in this study resulted in a comparative evolutionary dynamic analysis of the Cambodian RABVs including the estimation of the rate of nucleotide substitution per site, the TMRCA, and the MCC analysis of the complete coding region of Ngene sequences. Our analysis suggests that the SEA1 lineage had a most recent ancestor in common with the newly established SEA2 lineage. The two lineages apparently diverged around 1950. Interestingly, the SEA2 lineage shared a common ancestor with the Lao viruses around 1953, suggesting that either this lineage was first introduced to Laos from Cambodia and then was largely outcompeted in Cambodia by the SEA1 lineage, or that sporadic introductions of the SEA2 RABVs took place from Laos to Cambodia during the 1960s. It is possible that this lineage is more widespread in Cambodia than it appears, but additional studies are required to confirm this hypothesis. During the conflict period of the Khmer Rouge regime, stray dogs almost disappeared due to starvation and human consumption, which could also have impacted the circulation of the RABV lineages in the country. Based on the limited number of sequences available for this analysis, we established that viruses from Myanmar diverged from the other viruses of the major SEA clade approximately in 1942. Interestingly, the viruses of Chinese origin included in this study had a TMRCA around 1734 and since then seem to have evolved independently from the SEA clade. The present study identified three phylogenetically distinct lineages inside the major South-East Asian phylogenetic clade, with distributions representatives from the five countries in the region. Because the SEA1 lineage was the most widespread in Cambodia, the detection of a new SEA2 lineage in the country could have resulted either from a recent introduction of RABVs from abroad (e.g. Laos), or from the identification of a group of less prevalent viruses. Further molecular and epidemiological investigations are required to answer this question. To elucidate the distribution of the local lineages of the SEA clade and their roles in the molecular epidemiology of rabies, namely the dynamics and spread

Figure 1 .

 1 Figure 1. Geographical origins of the Cambodian rabies viruses sequenced for the purpose of the present study. The map was built with ArcMap 10 package (ESRI, Redlands, CA, USA).

Figure 2 .

 2 Figure 2. Evolutionary relationships of South-East Asian rabies viruses based on completeN-gene sequences. The tree was built using Maximum Likelihood method based on Tamura 3-parameter method[START_REF] Tamura | Estimation of the number of nucleotide substitutions when there are strong transition-transversion and G + C-content biases[END_REF] to compute the evolutionary distances. The tree with the highest log likelihood (-7633.4158) is shown. The percentage of trees in which the associated taxa clustered together is shown next to the branches. The tree is drawn to scale, with branch lengths measured in the number of substitutions per site. The analysis involved 192 nucleotide sequences. There were a total of 1346 positions in the final dataset. Bootstrap values lower than 70 are hidden.

Figure 3 .

 3 Figure 3. A. Maximum clade credibility tree of 190 complete coding region of nucleoprotein gene sequences of rabies viruses. The estimated TMRCA are indicated next to each node. Branches are color-coded by country of origin: dark blue for Cambodia, pink for Vietnam, red for Thailand, green for Laos, sky-blue for Myanmar and purple for China. The horizontal axis is shown as the number of years from the most recent sample. B. Population dynamics of rabies viruses (complete coding region of nucleoprotein gene sequences). The Bayesian skyline plot shows changes in relative genetic diversity (X axis) and time (Y axis) is shown as the number of years from the most recent sample. The black solid line represents the mean estimate of Neτ, while the 95% HPD intervals are shown in grey solid line. RABV sequence JN786877was excluded from this analysis (absence of information regarding the date of isolation).

Table 1 .

 1 PCR and sequencing primers used in this study

		Name of
	No	PCR primers Sequence 5' -3'	References
			PCR primers
	1	PV 05m	ATG ACA GAC AAY YTG AAC AA
	2	PV 09	TGA CCA TTC CAR CAR GTN G
	3	PV 08	GGT CTG ATC TRT CWG ARY AAT A	Dacheux L. et al., 2008
	4	β -Taq1	TCACCCACACTGTGCCCATCTACGA
	5	β -Taq2	CAGCGGAACCGCTCATTGCCAATGG
			Sequencing primers
	1.	RVN-71F	ATGGATGCCGACAAGATTGTATTC
	2. 3.	RVN-1118R RVN-1091F	GAATTCCTCTCCCAGATAGCC CTAGGGGGCTATCTGGGAGA	Yang et al., 2011
	4.	RVN-1562R	CGGCCAGACCGGCTCTAACAC

Acknowledgments

The authors would like to thank staff of the IPC Rabies Centre and International Centre for Vaccinations and the staff of the IPC Virology Unit.

6.

5.

APPENDIX

APPENDIX Supplementary Figure 1. Evolutionary relationships of rabies viruses from South-East Asia based on partial N-gene sequences. The tree was built using Minimum Evolution method [START_REF] Rzhetsky | A simple method for estimating and testing minimum evolution trees[END_REF]. The optimal tree with the sum of branch length = 304.71612549 is shown. The percentage of replicate trees in which the associated taxa clustered together in the bootstrap test (1000 replicates) are shown next to the branches [START_REF] Felsenstein | Confidence limits on phylogenies: An approach using the bootstrap[END_REF]. Boot strap values lower 70 are hidden. There were a total of 381 positions in the final dataset.

Supplementary Table 2.

Information on the reference sequences

No

Strain Name Country Year Species GenBank accession #