Sabine Grimm

PhD Carles Escera
email: cescera@ub.edu

Lavinia Slabu

Jordi Costa-Faidella

Electrophysiological evidence for the hierarchical organization of auditory change detection in the human brain

published or not. The documents may come

Electrophysiological evidence for the hierarchical organization of auditory change detection in the human brain

Sabine Grimm, Carles Escera, Lavinia Slabu, Jordi Costa-

The ability to detect new events in the acoustic environment is vitally important, as they might call for a prompt adaptive response. This requires that regularities in the acoustic input are modeled and kept in memory, so that deviant or contextually new stimuli violating the regularity representations can be detected. Processes of deviance detection have traditionally been associated with a particular component, the mismatch negativity (MMN; [START_REF] Escera | The mismatch negativity 30 years later: How far HaveWe come[END_REF][START_REF] Näätänen | Early selective-attention effect on evoked potential reinterpreted[END_REF][START_REF] Näätänen | The mismatch negativity. where is the big fish[END_REF] of the human event-related potential (ERP). This brain response is usually obtained with the auditory oddball paradigm comparing activity elicited by a frequently repeated stimulus (standard) to that elicited by an interspersed rare stimulus containing a feature variation (deviant). In that way, MMN can be obtained for violations of simple feature rules, as for example in the case of frequency, location, or intensity deviants, but it is also elicited for sounds violating more complex regularities (e.g. phonetic contrasts, abstract regularities defining the relationship between sounds etc.; for an overview, see [START_REF] Picton | Mismatch negativity: Different water in the same river[END_REF]. Thereby an MMN response is generated at 100 -250 ms from deviance onset by sources located bilaterally in the supratemporal brain region in the vicinity of auditory cortex [START_REF] Alho | Cerebral generators of mismatch negativity (MMN) and its magnetic counterpart (MMNm) elicited by sound changes[END_REF][START_REF] Näätänen | Generators of electrical and magnetic mismatch responses in humans[END_REF][START_REF] Maess | Localizing pre-attentive auditory memory-based comparison: Magnetic mismatch negativity to pitch change[END_REF]. Additional prefrontal contributions have been reported in several studies [START_REF] Deouell | The frontal generator of the mismatch negativity revisited[END_REF].

The underlying processes involve the modeling and storage of the acoustic regularities and can not merely be explained by different states of refractoriness of feature-specific neurons responding to the standard or deviant [START_REF] Näätänen | Memory-based or afferent processes in mismatch negativity (MMN): A review of the evidence[END_REF]; see also [START_REF] Jääskeläinen | Human posterior auditory cortex gates novel sounds to consciousness[END_REF]. This has been shown by means of a controlled oddball paradigm for a variety of deviant types including location [START_REF] Schröger | Mismatch response of the human brain to changes in sound location[END_REF], frequency [START_REF] Jacobsen | Is there pre-attentive memory-based comparison of pitch?[END_REF], and duration [START_REF] Jacobsen | Measuring duration mismatch negativity[END_REF] deviants.

In the controlled paradigm, the deviant stimulus from the oddball block is compared to a physically identical sound occurring with the same probability as the deviant in a context of different randomly presented equiprobable stimuli. Thus, the differential response is ensured not to be due to the differences in stimulus probability and associated differences in the state of refractoriness of neural populations, but is reflecting "true" deviance detection based on a regularity representation stored in auditory sensory memory.

Nevertheless, our understanding of the neural mechanisms underlying auditory deviance detection is still fragmentary. A step forward in this direction has been recently provided by studies of single-unit recordings in anesthetized animals. Indeed, the majority of neurons of the cat's primary auditory cortex exhibit a property termed stimulus-specific adaptation (SSA), that is, they reduce significantly their discharge rate after a few repetitions of the standard tone, but show fast robust responses to novel1 stimuli that slightly differ in their feature properties from the standard [START_REF] Ulanovsky | Processing of low-probability sounds by cortical neurons[END_REF][START_REF] Ulanovsky | Multiple time scales of adaptation in auditory cortex neurons[END_REF].

Despite sharing similar characteristics with the human MMN, the early latency of these noveltyspecific neural responses (circa 20 ms) suggests that they are not directly equivalent, but rather lie upstream of MMN generation (von der Behrens, Bäuerle, Kössl, & Gaese, 2009). Furthermore, very recent single-unit studies have shown that neurons in the inferior colliculus in the rat [START_REF] Malmierca | Stimulus-specific adaptation in the inferior colliculus of the anesthetized rat[END_REF][START_REF] Pérez-González | Novelty detector neurons in the mammalian auditory midbrain[END_REF] and in the barn owl [START_REF] Reches | Stimulus-specific adaptations in the gaze control system of the barn owl[END_REF] and neurons in the medial geniculate body of the thalamus [START_REF] Anderson | Stimulus-specific adaptation occurs in the auditory thalamus[END_REF][START_REF] Antunes | Is there stimulus-specific adaptation in the auditory thalamus?[END_REF] exhibit similar SSA to that found in the cat A1, suggesting that deviance detection in the auditory system can be found even before the information reaches the auditory cortex.

A comprehensive interpretation of these animal and human results is suggestive of two important aspects regarding the auditory system: 1) that deviance detection is a key principal expanding along the auditory pathway from the lower levels of the brainstem to high-order areas of the cerebral cortex; 2) that the generation of the MMN recorded from the human scalp is the consequence of a cascade of deviance detection processes at these different levels. Yet, a unified picture of these two lines of research, in humans and animals, is missing.

In the present study we aimed to test the hypothesis of a "pervasive auditory novelty system" by investigating processes of auditory change detection in humans on multiple time scales including the time range corresponding to the novelty responses observed in auditory cortex of animals.

Therefore, we measured aside from the MMN component an earlier portion of the ERP, the human Middle Latency Response (MLR) during a controlled oddball paradigm. The MLR is characterized by a sequence of waveforms in the range 12-50 ms from sound onset, labeled as P0, Na, Pa, and Nb [START_REF] Picton | Human auditory evoked potentials. I. evaluation of components[END_REF], composed of activation in subcortical, and primary and secondary auditory cortices [START_REF] Deiber | Sequential mapping favours the hypothesis of distinct generators for Na and Pa middle latency auditory evoked potentials[END_REF][START_REF] Liegeois-Chauvel | Evoked potentials recorded from the auditory cortex in man: Evaluation and topography of the middle latency components[END_REF][START_REF] Yvert | Multiple supratemporal sources of magnetic and electric auditory evoked middle latency components in humans[END_REF]Yvert et al., 2002). The MLR and MMN responses were measured to frequency deviants (800 Hz and 3730 Hz, in separate conditions) occurring in an oddball block and the responses were compared to those elicited by physically identical stimuli when they had the role of a standard in a "reversed" oddball block, and when occurring equiprobably amongst other four low-probability tones (control condition; Figure 1).

The setup of the EEG recordings was tailored to provide the possibility to extract the MLR components P0, Na, Pa, Nb and the long-latency component MMN in parallel analyses.

Our results revealed "true" deviance detection at a latency of 40 ms in humans, i.e., by the Nb waveform of the MLR, supporting the idea of a multistage comparison system for change detection along the auditory pathway in humans.

Methods

Participants

Twenty healthy, normal-hearing students (18-31 years, 11 female) participated in the experiment for payment (6 € per hour). All participants had normal hearing at both ears with a mean hearing threshold below 25 dB tested for the five frequencies used in the experiment. The experimental protocol was approved by the Ethical Committee of University of Barcelona, and was in accordance with The Code of Ethics of the World Medical Association (Declaration of Helsinki).

Participants gave written informed consent before the experiment.

Materials

Auditory sequences were composed of pure sinewave sounds of 50 ms duration including a 5 ms rise and a 5 ms fall time. Sounds were presented binaurally via headphones at an intensity level of 50 dB above the individual hearing threshold as measured at the beginning of each experiment. The stimulus onset-to-onset interval was set to 293 ms. The oddball sequence contained a frequent standard sound occurring with a probability of 0.80 and a rare frequency deviant occurring randomly with a probability of 0.20. In one of two conditions, the standard frequency was 1200 Hz and the deviant frequency was 800 Hz (low frequency condition). In the second, high frequency condition, the standard frequency was 2580 Hz and the deviant frequency 3730 Hz. Corresponding to each of the oddball blocks, a reversed sequence was introduced in which the roles of deviant and standard stimuli were switched. Additionally, a control condition was presented intermixing randomly tones of five different frequencies (800, 1200(800, , 1780(800, , 2580, 3730 Hz), 3730 Hz), each occurring with a probability of 0.2 from which the two extreme stimuli were taken as control tones for the low and high frequency condition, respectively (see Figure 1). This was done to preclude refractoriness confounds (i.e., the assumption that the deviant stimulus will elicit a stronger response per se, as it occurs with a much smaller probability than the standard and therefore would find feature-specific neural populations in a less refractory state than those responding to the standard stimulus [START_REF] Schröger | Mismatch response of the human brain to changes in sound location[END_REF].

Procedure

During the experiment, subjects were seated comfortably in an electrically shielded and soundattenuated chamber. They were instructed to relax and to watch a silent movie with subtitles ignoring the auditory stimulation. In total, 920 trials per deviant, standard and control stimulus were delivered. The three conditions (oddball, reversed oddball, and control) were split into a total number of 16 blocks of approximately 5 minutes each, which were presented in random order. After each block, subjects had a short break allowing for movements, after every fourth block a five minutes break for rest was introduced.

============================ FIGURE 1 ==============================

Electroencephalographic Recording

The electroencephalogram (EEG) was recorded continuously from 9 tin scalp electrodes referenced to an electrode placed on the tip of the nose. Electrodes were mounted according to the 10-20 system using an elastic cap (ECI Electro Cap) at the positions F3, F4, Fz, C3, C4, Cz, Pz and the left and right mastoids. Additionally, eye movements were measured bipolarly by two electrodes placed above and below the right eye (vertical electrooculogram, VEOG) and two electrodes placed at the outer canthi of each eye (horizontal electrooculogram, HEOG).

The electrode signals were amplified using a SynAmps amplifier (Neuroscan), online bandpassfiltered from 0.05 to 500 Hz and digitized with a sampling rate of 2500 Hz. Off-line data were re-referenced to the left mastoid.

EEG analysis

For the analysis in the long-latency range of the ERP, data was filtered off-line with a Kaiserwindowed sinc bandpass filter (beta=5.658) from 0.6 to 35 Hz. Epochs of 400 ms including a 100 ms pre-stimulus baseline were averaged separately for the deviant, standard and control stimulus in the two frequency conditions. Any trial with an amplitude variation larger than 80 µV was excluded from further analysis. Mean amplitudes of MMN were extracted at the electrode Cz from a 30-ms time window around the grand-average peak latency ranging from 90 to 120 ms as derived from the difference waveforms between deviant and standard, and deviant and control stimuli respectively. A repeated measures ANOVA including the factors Stimulus Type (deviant, standard, control) and Tone Frequency (low, high) was calculated on the mean amplitudes in the MMN time window. If appropriate, pairwise differences between single levels of factors were tested applying repeated measures contrasts.

For the MLR analysis, data was filtered with a Kaiser-windowed sinc bandpass filter (beta=5.658) from 15 to 200 Hz. Epochs of 150 ms including a 50 ms pre-stimulus baseline were averaged. Any trial with an amplitude variation larger than 80 µV was excluded from further analysis. The components P0, Na, Pa, and Nb of the MLR were extracted for the three stimulus types (deviant, standard, controls) and the two frequency conditions. Individual peak latencies were derived from the largest peak in the time windows 7-17 ms (P0), 19-29 ms (Na), 26-36 ms (Pa), and 37-47 ms (Nb), respectively. As mean amplitudes are known to be more reliable than peak amplitudes (Picton et al., 2000), analyses of the components' amplitudes were based on the mean voltage measured in a 4-ms time window centered on the respective mean grand-average peak latency elicited by the deviant, standard and control stimuli at the electrode Cz. For the high frequency stimuli, P0 peaked in the grand-average waveforms at 10 ms, Na at 23 ms, Pa at 28 ms, and Nb at 38 ms. For the low frequency stimuli, the respective grand-average peak latencies were at 12 ms (P0), 24 ms (Na), 32 ms (Pa), and 42 ms (Nb) after tone onset.

For each component a repeated measures ANOVA including the factors Stimulus Type (deviant, standard, control) and Tone Frequency (low, high) was calculated on MLR peak latencies and mean amplitudes. If appropriate, pairwise differences between single levels of Stimulus Type were tested applying repeated measures contrasts.

A result was considered significant when p < 0.05 using a two-tailed analysis. Bonferroni correction was used for all multiple pairwise contrasts as well as for testing on multiple MLR components. The partial eta squared (η p 2 ) is reported in addition to F-and p-values.

Results

ERPs recorded in 20 participants displayed the typical series of components in response to the standard, deviant, and control tones in the low (Figure 2 a) and in the high (Figure 2 b) frequency conditions. Focusing on the long-latency range of the ERPs, a prominent MMN was obtained peaking at about 105 ms after tone onset. In the MLR, the characteristic P0-Na-Pa-Nb complex was displayed in all experimental conditions (Figure 2, bottom). MLR peak latencies are given in Table 1. Table 2 shows the mean amplitude values for each MLR component in the 4-ms latency windows around the grand-average peaks as given above.

Long-latency components of the ERP

First, data were analyzed with respect to the long-latency range of ERP contrasting evoked responses to standard, deviant, and control stimuli for both the low and high frequency tones. In the time window around the MMN peak from 90 to 120 ms, deviant ERPs displayed a sharp negative potential, this being of less negative amplitude for control ERPs and even positive for standard ones. Amplitude differences for the three stimulus types were statistically significant in 1.35, corrected p = 1.000, η p 2 = 0.07; cf. Figure 3). Post-hoc repeated measures contrasts confirmed differences in the Nb latency window between deviant and standard responses (F(1,19) = 9.98, corrected p = 0.010, η p 2 = 0.34) and between deviant and control responses (F(1,19) = 6.57, corrected p = 0.038, η p 2 = 0.26). Amplitude modulations through Stimulus Type did not interact with the effects of Tone Frequency for any of the components (F's < 1.1, corrected p's = 1.000, η p 2 < 0.06).

Discussion

The results of the present study have revealed that "true" auditory deviance detection can take place in humans as early as 40 ms after the onset of a deviant feature presented in an otherwise repetitive sequence of standard stimuli. This was supported by the fact that a characteristic amplitude modulation of the MLR component Nb (peaking at about 38-42 ms) was obtained depending on the deviant status of a stimulus. The Nb response was larger for a stimulus when it occurred as a frequency deviant than the response elicited by the same stimulus when it occurred in the role of a standard or a control tone. In the later portion of the ERP, we additionally obtained a clear MMN peaking at about 105 ms after sound onset when contrasting deviant against standard and control responses.

The MMN is the typical marker of deviance detection known to reflect a memory-based process of comparing incoming stimuli with an internal model derived from the regularities in the previous stimulation. A respective control condition is required [START_REF] Schröger | Mismatch response of the human brain to changes in sound location[END_REF] to ensure that any deviance-related modulation obtained in an oddball paradigm is truly memorybased and not merely reflecting differences in the response strength of feature-specific neural populations (that are more refractory in case of a repeatedly presented standard frequency than in case of a rarely presented deviant one). Of similar importance is the control for physical stimulus properties. Particularly the MLR is a sequence of components whose latencies and amplitudes are systematically influenced by physical stimulus characteristics [START_REF] Picton | Human auditory evoked potentials. I. evaluation of components[END_REF]; as reflected in the latency and amplitude differences for the MLR components Pa, Na and Nb between low and high frequency tones in the present study, which is congruent with previous reports in the literature [START_REF] Kraus | Color imaging of the human middle latency response[END_REF][START_REF] Thornton | Effects of halothane or enflurane with controlled ventilation on auditory evoked potentials[END_REF].

Taking into account these points, the MMN and the enhancement of Nb obtained in the present study are due to the informational status carried by a stimulus in its respective context and thus reflect true deviance detection rather than confounds by stimulus properties or refractoriness.

Thus, it can be concluded that besides the typical MMN component also the middle-latency portion of the ERP (which has so far barely been analyzed in studies applying the oddball paradigm) is sensitive to stimulus novelty. This goes in line with a few recent results challenging the long-held belief that early auditory processing reflected in the MLR solely depends on the physical properties of an incoming stimulus. Sonnadara and colleagues (2006) reported an enhanced Na component of the MLR peaking at 25 ms after sound onset for location deviants compared to standard click sounds thus indicating an early effect of stimulus rareness. Yet, the lack of a respective control condition does not ascertain conclusively whether those results can be attributed to "true" deviance detection. Furthermore, modulations of MLR amplitudes by sensory gating [START_REF] Müller | Suppression of the auditory middle-latency response and evoked gamma-band response in a paired-click paradigm[END_REF], self-initiation of a stimulus [START_REF] Baess | Attenuated human auditory middle latency response and evoked 40-hz response to self-initiated sounds[END_REF], task requirements [START_REF] Woldorff | Modulation of early auditory processing during selective listening to rapidly presented tones[END_REF] and even by sound segregation processes [START_REF] Dyson | Representation of concurrent acoustic objects in primary auditory cortex[END_REF]) have been shown, altogether underlining the complex nature of auditory processing already in its initial phase.

From animal studies we have indication that deviance detection can arise already at this early steps of auditory processing. For the first time, [START_REF] Ulanovsky | Processing of low-probability sounds by cortical neurons[END_REF][START_REF] Ulanovsky | Multiple time scales of adaptation in auditory cortex neurons[END_REF] comprehensively described the activity of novelty neurons in the primary auditory cortex (A1) of the cat. These neural responses share a variety of properties with the MMN for which they have been regarded its single-neuron correlate. Both the firing of novelty neurons and MMN are pre-attentive responses whose magnitude is inversely related with the deviant probability and positively related with the degree of deviance, both show already local sequence effects, their latencies are similarly influenced by deviant probability and both responses are localized to the auditory cortex [START_REF] Nelken | Mismatch negativity and stimulus-specific adaptation in animal models[END_REF].

On the other hand, however, there are remarkable differences in timing between the firing onset of novelty units (at about 20 ms from stimulus onset; [START_REF] Pérez-González | Novelty detector neurons in the mammalian auditory midbrain[END_REF][START_REF] Ulanovsky | Multiple time scales of adaptation in auditory cortex neurons[END_REF] and the peak latency of the MMN which is contradicting the view that the first directly accounts for the latter [START_REF] Von Der Behrens | Correlating stimulus-specific adaptation of cortical neurons and local field potentials in the awake rat[END_REF]. Therefore the activity of novelty neurons has been interpreted as a change detection process in PAC that lies upstream of later MMN generation. The modulation of Nb by stimulus deviance found here, exemplifies that those "upstream" activities also exist in the human auditory system and that they can be identified in the ERP with a respectively tailored setup. It remains open whether the deviance-related Nb enhancement is more directly linked to the activity of novelty neurons in PAC. At least the two share partly similar origins as the transition of components Pa to Nb is supposed to be generated by cortico-cortical connections mediating auditory information from PAC to the superior temporal gyrus (STG; Yvert et al., 2002). Yet, only the simultaneous use of the different techniques in future studies might permit to disambiguate the temporal relationship between single-cell firing and the potentials measured over larger auditory fields.

Eventually, the present study confirms that deviance detection is implemented on multiple levels during auditory processing in the human brain. In order to integrate these results into a conceptual framework, we here propose that deviance detection and the underlying processes of modeling invariant input are a pervasive property of the auditory system, expanding from lower levels along the auditory pathway to high-order areas of the cerebral cortex. This property allows us to react quickly to new events in our environment and crucially shapes our perception by sharpening its sensitivity to changes in the incoming information. Within this framework, it can be assumed that the generation of MMN is a consequence of a cascade of deviance detection processes occurring at hierarchically lower levels. Accordingly, we can hypothesize that this function is organized in a hierarchical manner, so that deviance-related responses to simplefeature changes as used here are detected at the lower levels of the novelty system's hierarchy whereas more complex levels of regularity will be encoded in higher levels and thus in the latency range of MMN only. Alternatively, one could speculate about more profound differences in the functional significance of the two levels of deviance processing, with the earlier possibly being related to a mechanism of auditory predictions (compare, [START_REF] Bendixen | I Heard That Coming: Event-Related Potential Evidence for Stimulus-Driven Prediction in the Auditory System[END_REF][START_REF] Winkler | Modeling the auditory scene: predictive regularity representations and perceptual objects[END_REF], and the later level being related to regularity updating which has been proposed previously as one of the functional roles of MMN [START_REF] Winkler | Adaptive modeling of the unattended acoustic environment reflected in the mismatch negativity event-related potential[END_REF][START_REF] Winkler | Interpreting the mismatch negativity (MMN)[END_REF]. Without question, future studies are needed to clarify the functional significance of the different levels of deviance processing.

To summarize, we report here a new electrophysiological marker of auditory deviance detection that indicates the auditory system's rapid response to rare and unexpected sounds. This response resembles in origin and latency very much the recently described responses of so-called novelty neurons found in the cat's auditory cortex and complements our picture of the functional organization of the auditory system. The presence of different markers of deviant processing in two time ranges of the ERP in the present study strongly supports the idea of a hierarchically organized system serving auditory deviance detection.

Figure Captions Figure 1 .

 Captions1 Figure Captions

Figure 2 .

 2 Figure 2. Grand-average evoked potentials of 20 subjects (a) for the low frequency tones and (b)

Figure 3 .

 3 Figure 3. Bar charts showing mean amplitudes of the MLR components (left side) and MMN

  Acknowledgement: This work was supported by the program Consolider-Ingenio 2010 (grant number CDS2007-00012), the National Program for Fundamental Research (grant number PSI2009-08063), a Juan de la Cierva grant to Sabine Grimm (grant number JCI-2009-04401) by the Spanish Ministry of Science and Innovation, and the ERANET NEURON project PANS (grant number EUI2009-04086). The authors thank Marc Recasens for his help during data acquisition.Yvert, B., Fischer, C., Guenot, M., Krolak-Salmon, P., Isnard, J., & Pernier, J. (2002). Simultaneous intracerebral EEG recordings of early auditory thalamic and cortical activity in human. The European Journal ofNeuroscience, 16(6), 1146Neuroscience, 16(6), -1150. .

Table 1 .

 1 Mean of individual peak latencies for the MLR components P0, Na, Pa, and Nb. Standard Errors of Mean are given in parentheses.

		P0	Na	Pa	Nb
		Mean peak	Mean peak	Mean peak	Mean peak
		latency (SEM)	latency (SEM)	latency (SEM)	latency (SEM)
	Freq 800 Hz
	Sta	13.35 (0.49)	24.05 (0.44)	32.50 (0.48)	43.35 (0.45)
	Dev	13.05 (0.33)	24.50 (0.29)	33.00 (0.36)	43.20 (0.58)
	Con	12.45 (0.49)	24.65 (0.40)	31.70 (0.48)	43.60 (0.55)
	Freq 3730 Hz
	Sta	11.90 (0.66)	23.30 (0.51)	30.50 (0.69)	42.10 (0.70)
	Dev	13.15 (0.60)	22.60 (0.47)	29.35 (0.44)	40.85 (0.65)
	Con	12.05 (0.59)	22.85 (0.51)	30.50 (0.59)	40.80 (0.67)

Table 2 .

 2 Mean amplitudes of the MLR components P0, Na, Pa, Nb derived from a 4-ms latency windows centered around the peak latency in the grand-average waveforms. Standard Errors of Mean are given in parentheses.

	P0		Na	Pa	Nb
	Mean window amplitude	Mean window amplitude	Mean window amplitude	Mean window amplitude
	(SEM)	(SEM)	(SEM)	(SEM)
	Freq 800 Hz
	Sta	0.122 (0.042)	-0.198 (0.038)	0.237 (0.045)	-0.203 (0.055)
	Dev	0.182 (0.040)	-0.178 (0.038)	0.292 (0.061)	-0.285 (0.051)
	Con	0.171 (0.047)	-0.195 (0.043)	0.233 (0.062)	-0.253 (0.062)
	Freq 3730 Hz
	Sta	0.027 (0.044)	-0.079 (0.037)	0.067 (0.043)	-0.083 (0.046)
	Dev	0.125 (0.024)	-0.069 (0.033)	0.101 (0.039)	-0.222 (0.046)
	Con	0.082 (0.036)	-0.143 (0.029)	0.042 (0.032)	-0.149 (0.036)

Note, that the term novel is used here in a more general meaning referring to contextually new or deviant information.


Mihail D Croitoru
email: mihail.croitoru@u-bordeaux1.fr

A A Shanenko

C C Kaun

F M Peeters

Ultra-small metallic grains: effect of statistical fluctuations of the chemical potential on superconducting correlations and vice versa

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

Introduction

In the past two decades quasi-low-dimensional superconductor systems have become a very active field of research due to their rich variety of physical properties and potential applications, which are superior to their three-dimensional homogeneous counterparts [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15]. The technological advances in fabricating and controlling material structures of nanoscale size provided a whole range of opportunities for exploring experimentally these physical properties, induced by quantum confinement [16][17][18][19][20][21][22][23][24][25][26][27]. Quantum confinement results in the discretization of the electronic energy bands and its competition with the pairing correlations produces a nonuniform spatial distribution of the superconducting order parameter.

The discreteness of the single-particle energy levels due to electron confinement plays a fundamental role in the description of the superconducting state in a metallic grain. A reduction of grain size leads to an increase of the spacing between the levels. Back in 1959, Anderson argued that when the average level spacing δ, which is inversely proportional to the density of states at the Fermi level and the grain volume, becomes of the order of the BCS bulk gap 0 then superconductivity should disappear [28]. Since the celebrated experiments by Ralph, Black, and Tinkham [16,17] on Al superconducting grains, several groups have reported the fabrication of different metallic superconducting nanoparticles, some of them with crystalline perfection [18,21,25,26]. Experiments have demonstrated clear signatures of superconducting pairing up to δ ∼ 0 , whereas no such correlations survived for δ 0 . This stimulated several theoretical investigations on the critical level spacing in nanograins at which superconducting correlations cease [1,2,29].

A mean-field model of a superconducting grain with equidistant energy levels was considered by von Delft et al [1]. This model takes into account number-parity effects, which are of particular importance in small grains. They found that the superconducting gap vanishes abruptly at a certain critical value of the level spacing δ c , which depends on the parity-number of electrons in the grain, namely, δ o c = 0.89 0 in grains with an odd number of confined electrons, and δ e c = 3.56 0 in grains with an even number of electrons. Smith and Ambegaokar [2] extended the Delft's mean-field treatment of superconductivity in such small grains by taking into account the Wigner-Dyson distribution of energy levels as obtained within random matrix theory [START_REF] Mehta | Random Matrices[END_REF]. They found that for both odd and even cases randomness enhances pairing and that superconductivity persists beyond the critical level spacing as deduced from the equal level spacing model. It turns out that the critical level spacings at which pairing vanishes abruptly are increased to δ o c = 1.80 0 for an odd number of electrons, and δ e c = 13.81 0 for an even number of electrons. The treatments beyond the mean-field approximation [5][6][7] have shown that superconducting pairing persists beyond the critical level spacing, but only in the form of fluctuations, and the crossover between regimes of conventional and fluctuation-dominated superconductivity is smoothed out [5,6].

In both mean-field models the chemical potential was chosen to lie halfway between the last filled and first empty levels for an even number of electrons (the even case), and on the half-filled level for an odd number of electrons (the odd case), which is a typical situation in a nanoparticle, when electron-electron correlation is neglected. However, the position of the chemical potential is of importance because the density of states in the vicinity of µ strongly influences the superconducting correlations, and hence the final conclusion on the critical level spacing. The position of the chemical potential in its turn is affected by the correlations. To our knowledge, a thorough investigation of the influence of the superconducting pair correlations in a nanograin on the position of the chemical potential has not been reported in the literature.

Since in the most experiments only the average behavior of a large number of grains is obtained, the aim of the present paper is to study a statistical distribution of the chemical potential around its position, which it has in the absence of electron-electron correlations (henceforth, µ 0 -position).

As an example, we extend the mean-field model, developed in [2], by weakening the µ 0 condition, and study how this effect can ostensibly modify the critical level spacing.

The layout of our paper is as follows. In section 2, we investigate the statistical distribution of the chemical potential induced by the pair correlations in a superconducting nanograin. In section 3, we show how the grain critical level spacing is affected by the statistical distribution of µ. Finally, a short summary is given in section 4.

Grain size dependence of the chemical potential distribution

We describe the grain by the following reduced BCS Hamiltonian, where only the time-reversal states are paired

H = p,σ ε p a † p,σ a p,σ - p,q
V p,q a † p↑ a † p↓ a q↓ a q↑ .

Here p, q are integers numbering the single-particle energy levels ε p and the operator a p,σ (a † p,σ ) annihilates (creates) an electron in state p with spin σ . The interaction matrix element V qp is given by

V qp = g d 3 r |ϕ q (r)| 2 |ϕ p (r)| 2 , (2)
with g denoting the coupling constant and ϕ q (r) the single-electron wavefunction. The first term in equation ( 1) contains the single-electron energies, and the second term is the attractive (when g > 0) pairing interaction due to the exchange of virtual phonons. We have assumed that the electron-electron interaction is unaffected by quantum confinement and it is the same as in the bulk. In the bulk the real inter-electron potential is well approximated by a δ-function pseudopotential. Employing such a simplified interaction requires a regularization, which makes the matrix elements non-zero only between states within the Debye window around the Fermi surface.

For a given value of the chemical potential the parity projected BCS self-consistency equation is [15,[START_REF] Balian | [END_REF]

η q = |ε p -µ|<hω D V qp η p 2E p t p 1 - ηR π t -2 q (1 -t -2 p ) 1 + ηR π t -2 q , (3)
where t p = tanh(βE p /2) with E p the quasi-particle energy, which is given by a complicated expression but is well approximated by

E p = [(ε p -µ) 2 + | p | 2 ] 1/2 ,
R π /2 ≡ p t p and η is the parity index (η = -1, 1 for the odd and even cases, respectively).

For a given number N e of electrons in the grain the chemical potential µ is determined from equation

N e = p f p 1 - ξ p E p + 1 -f p 1 + ξ p E p , (4)
with ξ p = (ε p -µ) and f p = (1 -t p )/2. The coupled set of equations ( 3) and ( 4) were solved numerically by iterations until full self-consistency was reached. At the first iteration the value of the chemical potential in the absence of superconducting correlations is taken. Numerical calculations were performed with parameters typically for Sn, namely, hω D /k B = 195 K and gN(0) = 0.25, with N(0) the bulk density of states at the Fermi level, and we used the bulk electron density n e = 148 nm -3 . In the rest of the paper we consider a nearly cubic grain with L x = D + a + ν 1 , L y = Da + ν 2 and L z = D + b, while varying the thickness of the grain D. Here a = 0.2 nm; ν 1 , ν 2 are the random deviates with a normal (Gaussian) distribution used to generate a large sample of single-electron level configurations in a grain. The standard deviation for the normal distribution was chosen to be 0.1 nm. b is fixed to a value that keeps the volume of the grain constant for given D constant, i.e. b = (D + a)(Da)D/L x L y . Parallelepiped-shaped grains were taken in order to avoid strong degeneracy of the single-electron levels typical for the highly symmetric samples. The single-electron wavefunction in such a grain reads

ϕ p (x, y, z) = l=x,y,z 2 L l sin πj l r l L l , (5)
where p = {j x j y j z }. In this case the interaction matrix element, equation ( 2), acquires the form

V qp = g l=x,y,z v q l ,p l , (6)
where

v q l ,p l =        3 2L l , if q l = p l , 1 L l , if q l = p l . (7)
The odd-N BCS solution with lowest energy assigns a specific role to the twofold-degenerate individual level whose energy ε p is closest to the chemical potential µ. Throughout this paper, any quantity pertaining to this level is labeled by the index p = 0, and the level with second smallest quasi-particle energy is denoted by the index p = 1. For zero temperature, q = 0 and η = -1, equation ( 3) is simplified to

q =0 = p(p =0) V qp p 2E p , (8)
and (for q = 0)

q=0 = p(p =0,1) V 0,p p 2E p . (9)
For convenience, we no longer mention explicitly that only single-electron states with |ε p -µ| < hω D are taken into account; i.e., the odd-ground state is obtained by creating a p = 0 quasi-particle above the BCS vacuum. The only formal difference between equation ( 8) and the usual BCS model is the suppression of the term with p = 0. Nevertheless, this difference may lead to a significant decrease of the pairing [1,[START_REF] Balian | [END_REF],

because the p = 0 quasi-particle prohibits Cooper pair scattering involving p = 0. In nuclear physics, the above approximation is usually referred to as the blocking (Soloviev) approximation [32]. The gap q=0 does not enter in equation (8), which determines all the other gaps at sufficiently low temperatures. However, the presence of such a gap associated with the blocked level allows one to relax the requirement of µ = µ 0 = ε 0 for the odd case [33]. We can expect that µ-ε 0 can be of the order of q=0 ∼ 0 . Indeed, this is seen from our numerical simulations.

Figures 1(a) and (b) show the probability distribution of the chemical potential in odd grains with D = 2.107 nm without and with the blocking effect, respectively. The numerical simulations were performed with N = 2.0 × 10 4 energy level configurations and with energy steps E = 0.025 meV (the number of discrete nodes for sampling the energy space is 200). The distribution is scaled by its maximum value. In the absence of superconducting correlations in the grains all single-electron level configurations result in a chemical potential pinned to the half-filled level µ = ε 0 for an odd number of electrons. The zero energy in the figure is the position of this level, µ 0 = 0. One can see from figure 1(a) that switching on the inter-electron correlations leads to a chemical potential probability distribution that is different from the δ(µµ 0 )-function. However, as seen from figure 1 to a decrease of the relative height of the central peak in the probability distribution. For large grains the redistribution of the chemical potential due to the pair correlation is negligible and when calculating superconducting characteristics for the chemical potential one can use its value for the grain without correlation as an acceptable approximation. The reason is that in this case there are a lot of single-electron levels in the Debye window, which allows the system to conserve the particle-hole symmetry without a large perturbation of the chemical potential position.

Figures 3(a)-(c) illustrate the probability distribution of the chemical potential in the even grain case with thicknesses D = 2.33, 3.448, and 8.0 nm. The behavior of the probability distribution is the same as in the odd case. The only noticeable difference is that the height of the redistribution does not approach zero when reducing the width of the sample (at least for the given range of the grain thicknesses).

To investigate the size dependence of the probability redistribution more qualitatively we fit part of the distribution, from which the central peak is removed, by a Gaussian function of the form The dependence of the standard deviation of the Gauss distribution of the part of the original distribution is strongly correlated with the strength of the pair correlation, which is described by the mean value of the superconducting order parameter. This is seen from figure 5, where the order parameter as a function of the grain thickness is shown [15]. To plot this figure we performed calculations for a parallelepiped-shaped grain with L x = D + a, L y = Da and L z = D, while varying the thickness of the grain D. Comparing this figure with the lower panel of figure 4 we can see that for D < 6-8 nm the average of the order parameter is approximately equals to 2σ .

G(µ, σ ) = B σ π 2 e µ-µ 0 σ 2 . ( 10

Critical level spacing

Let us now investigate the possible consequences of the chemical potential distribution on the critical level spacing in the superconducting grains.

The odd case (η = -1): for the purpose of simplicity in the rest of our work we assume that quantum confinement does not modify the matrix elements V qp , i.e. we disregard the spatial variations of the order parameter [9,10], and V qp 's are taken the same as in the bulk, namely V qp = g V = λδ, with λ the dimensionless coupling. This assumption results in p (δ, T) = (δ, T). Now, keeping this in mind and using equation ( 8), for k B T/ 0 1 we obtain

1 = λ δ p(p =0) t p 2E p . ( 11
)
Adding and subtracting the energy ε 0 of the single-electron level p = 0, which is closest to the chemical potential µ, we obtain 1 = λδ

2 hω D -hω D dω t ω p =0 δ ω -ε p + ε 0 (ω + ε 0 -µ) 2 + | (δ, T)| 2 . ( 12
)
The position of the electron energy levels in a grain strongly depends on the boundary conditions and the interaction of the electrons inside a grain. However, disorder related to irregularities of the shape of individual grains, or internal defects inside a grain, plus the interaction with charged impurities in the insulating substrate, make a calculation of the energy levels a daunting task. To overcome this difficulty one usually performs a statistical description of the grain characteristics, instead of calculating the precise position of the energy levels in a specific dot. In this case the energy spectrum is characterized by the mean level spacing and by the level fluctuations [START_REF] Mehta | Random Matrices[END_REF]34]. Thus we average equation ( 12) over disorder and shape irregularities

1 = λ 2 hω D -hω D dω t ω R ω π δ (ω + ε 0 -µ) 2 + | (δ, T)| 2 , (13)
where we introduce

R(ε)/δ = p δ ε -(ε p -ε 0 ) . ( 14
)
As shown in [2], this average can be related to the two-level correlation function (TLCF), describing the level fluctuations in the system. Since we consider the system without magnetic or spin-orbit interactions, which possess time-reversal and spin-rotation invariance, the TLCF given by RMT is that for a Gaussian orthogonal ensemble [START_REF] Mehta | Random Matrices[END_REF]

R(x) = 1 - sin 2 x x 2 - d dx sin x x ∞ x dt sin t t . (15)
Inserting equation ( 15) into (13), we arrive at the basic equation of the Smith-Ambegaokar model [2].

In the next step, we treat the position of the chemical potential statistically as was discussed in section 2 4 . As follows from that discussion, we have to model the statistical distribution of the chemical potential as the sum of two functions: one describing the central peak and another one for the redistributed part. However, since we are only interested in qualitative effects due to the redistribution of the chemical potential, we will neglect the central peak and approximate the redistribution part by the Gaussian function. Our approximation will result into an upper limit for the critical level spacing. The opposite situation, when we take into account only the central peak, leads to the Smith-Ambegaokar model.

So, the chemical potential is taken normally distributed around ε 0 . In addition, we assume that (ε p -ε 0 ) and (ε 0 -µ) are independent random variables. Then, we obtain

1 = λ hω D 0 dω R πω δ × ∞ -∞ dy t ω G (y, σ ) (ω + y) 2 + | (δ, T)| 2 , (16)
where G (y, σ ) is the Gaussian distribution function of the chemical potential given by equation (10).

The critical level spacing δ o c for the odd case is defined as the solution of equation ( 16) in the limit (δ, T) → 0, thus

ln 2 hω D 0 = π hω D /δ o c 0 dx R (x) ∞ -∞ dy t x,y G (y, σ ) x + πy/δ o c , (17)
where we introduce the new variable x = πω/δ and take into account that λ -1 = ln(2 hω D / 0 ) with 0 = (0, 0). Here t x,y = tanh(

= π hω D / δ o c 0 dx R(x) x . (18)
This can be solved analytically resulting in the odd-case result of Smith and Ambegaokar [2], i.e., δ o c = πe γ +π 2 /16-7/4 0 ≈ 1.8 0 , with γ the Euler constant.

The integrations in equation ( 17) are done numerically, from which we obtain the critical level spacing as δ o c = δ o c (σ, T). Figure 6 illustrates how the critical level spacing δ o c depends on the standard deviation for the distribution of the chemical potential. The mean value of the chemical potential in the odd case is the energy of the half-filled level ε 0 . Since, 4 The physics behind this is as follows. The distribution of energy levels is a function of a certain number of statistical variables. Since the quasiparticle energy is defined according to the expression

E p = [(ε p -µ) 2 + | p | 2
] 1/2 , the distribution of the chemical potential µ is a function of the same statistical variables (levels) plus an additional variable due to the order parameter. This additional variable is a function of the statistics of the system wavefunctions. In the GOA ensemble the statistics of the levels is independent of the statistics of the wavefunctions [35]. as it was shown in section 2, the deviation of the chemical potential from ε 0 is due to electron correlations (for T → 0) characterized by , the standard deviation in this figure is normalized by 0 . Note that fluctuations in the position of the chemical potential significantly enhance the critical level spacing, i.e., the regime of conventional superconductivity can persist in samples with smaller sizes. This effect can be explained as follows. The pair correlations are different from zero within a finite energy range around the chemical potential µ, and are strongest exactly at µ [START_REF] Balian | [END_REF]. However, the unpaired particle which occupies the level ε 0 forbids the formation of the pair state ( 0, 0), with 0 the time-reversed counterpart. Consequently, when the chemical potential is pinned to the level ε 0 we obtain the strongest decrease of the pairing correlations in the system. The redistribution of the chemical potential frees up phase space for the strongest pair scattering and so achieves a gain in interaction energy, and hence in the critical level spacing. This effect is suppressed by temperature because temperature diminishes the effect of Soloviev blocking. At finite-T the peak in the pairing correlations at the chemical potential level diminishes and, hence, the effect of the redistribution of the chemical potential is weakened for increasing T.

The even case (η = +1): the chemical potential is now distributed around the energy being halfway between the last filled level ε 0 and the first empty level ε 1 . The equation for the critical level spacing can be obtained as follows. We start with the chemical potential residing exactly at µ 0 = (ε 0 + ε 1 )/2, which yields

1 = λδ p t p 2 ε p -µ 0 + µ 0 -µ 2 + | (δ, T)| 2 . (19)
Then, for ε p = ε 0 and ε p = ε 1 we can write ε p -µ 0 = (ε p -ε 1 ) + (ε 1 -ε 0 )/2. Following similar arguments as for the odd case we found

1 λ = ω D 0 dω ∞ -∞ dz t ω R(ω π δ ) G (z, σ ) (ω + ε 1 -ε 0 2 + z) 2 + | (δ, T)| 2 + ∞ -∞ dz t ω G (z, σ ) ( ε 1 -ε 0 2 + z) 2 + | (δ, T)| 2 + ∞ -∞ dz t ω G (z, σ ) ( ε 0 -ε 1 2 + z) 2 + | (δ, T)| 2 .
The distribution of (ε 1 -ε 0 ) is given by the nearest level spacing distribution, P (y). This distribution is well approximated by the 'Wigner surmise' expression [START_REF] Mehta | Random Matrices[END_REF] P

(y) = y 2π exp - y 2 4π , (20)
where y = π (ε 1 -ε 0 ) /d. Assuming now that the distributions of ε p -ε 1 , (ε 1 -ε 0 ) and ε 0 +ε 1 2 -µ can be treated independently, we obtain the following equation for the critical level spacing δ e c in the even case: 1 λ = 2π This equation is solved numerically and the results for δ e c are shown in figure 7 as function of σ/ 0 . Now the mean value of the chemical potential is µ 0 . Notice that deviations of µ from µ 0 result in a monotonic increase of the critical level spacing δ e c with increasing standard deviation σ . This can be explained as follows. The pair correlations are strongest exactly at µ.

The redistribution of the chemical potential shifts one of the neighboring levels (either ε 0 or ε 1 ) towards the maximum of the pairing correlation function-or in other words shifts one of them towards the chemical potential-where the pair scattering is stronger. At the same time the pair scattering that involves another neighboring level becomes weaker. However, an increase of the correlations on the former level exceeds the decrease on the latter one, resulting in a total gain in interaction energy and, hence, in the critical level spacing. As expected, temperature diminishes this effect. Notice that for T = 0 and σ = 0 we obtain δ e c / 0 = 18.7, which is substantially larger than the approximate estimate δ e c / 0 ≈ 13.8 found by Smith and Ambegaokar [2].

Conclusion

Concluding, we have studied the influence of the superconducting pair correlations in a grain on the probability distribution of the chemical potential. It is shown that pair correlations result in a redistribution of the chemical potential, which depends on the grain size. In the odd case, for ultra-small grains the effect of the redistribution is sufficiently suppressed due to the blocking effect.

We also investigated the possible effect of the statistical fluctuations of the chemical potential on the superconducting correlations in nanoparticles. A statistical description of the chemical potential is justified by the presence of a pairing gap in superconducting grains. The deviation of the chemical potential from a fixed value typical for normal metals frees up the phase space that becomes available for the scattering of Cooper pairs. As a result, the critical level spacing increases.

Figure 1 .

 1 Figure 1. The chemical potential probability distribution for a superconducting Sn ensemble of grains with thickness D = 2.107 nm: (a) without Soloviev blocking; (b) with Soloviev blocking.

  (b), the blocking of the p = 0 level weakens the effect of the redistribution of the chemical potential around µ 0 by weakening the pair correlations. It has a very sharp peak in the very vicinity of its average value (for |E| = |µ -µ 0 | < E). Nevertheless, for this particular case of D = 2.107, 66.5% of the configurations resulted in a chemical potential in the range |µ -µ 0 | > E. Figures 2(a)-(c) illustrate the probability distribution of the chemical potential in odd grains with thicknesses D = 3.068, 5.190, and 7.5 nm. The inset in the upper panel is a zoom of the lower part of the probability distribution for D = 3.068 nm. It is seen that an increase of the grain size leads to a decrease of the distribution width and at the same time

Figure 2 .

 2 Figure 2. The chemical potential probability distribution for a superconducting Sn ensemble of odd grains with thicknesses (a) D = 3.068 nm; (b) D = 5.190 nm; (c) D = 7.5 nm.

Figure 3 .

 3 Figure 3. The chemical potential probability distribution for a superconducting Sn ensemble of even grains with thickness (a) D = 2.330 nm; (b) D = 3.448 nm; (c) D = 8.0 nm.

)π 2 is

 2 Then the evolution of the chemical potential distribution with grain thickness can be represented with the found dependence of the fitting parameters σ and B on the grain thickness. The result is shown in figures 4(a) and (b). The upper panel here illustrates the relative height of the Gauss distribution (the height A = B σ √ normalized by the height of the central peak in the original distribution). The curve is obtained by fitting our numerical results with the function y = a + bx c , with a = 0.0, b = 0.005 and c = 2.3 in the odd case and with a = 0.28, b = 6.6 × 10 -4 and c = 2.8 in the even case. The lower panel shows the standard deviation σ as a function of the grain size. The results were fitted by the curves with: a = 0.0, b = 9 and c = -2 (odd case), a = 0.0, b = 4.9 and c = -1.8 (even case).

Figure 4 .

 4 Figure 4. The relative height (panel (a)) and the standard deviation (panel (b)) of the normal distribution as a function of the grain size D.

Figure 5 .

 5 Figure 5. The mean value of the superconducting order parameter versus the grain thickness D.

  δ o c |). Equation (17) can only be solved numerically. When neglecting fluctuations of the chemical potential, i.e., G(y, σ = 0) = δ(y), this equation reduces to (at T = 0) ln 2 hω D 0

Figure 6 .

 6 Figure 6. The odd case: the critical level spacing δ o c versus the standard deviation σ of the chemical potential for different temperatures ranging from k B T/ 0 = 0 to 0.15.

Figure 7 .

 7 Figure 7. The even case: the critical level spacing δ e c versus σ for different temperatures from k B T/ 0 = 0 to 0.15.

1 x

 1 + y/2 -πz/δ e c .

Acknowledgments

This work was supported by the European Community under the Marie Curie IEF Action (Grant Agreement No. PIEF-GA-2009-235486-ScQSR), the Flemish Science Foundation (FWO-Vl), the Belgian Science Policy (IAP), and the ESF network INSTANS. MDC and AAS are grateful to A Vagov for stimulating discussions.


M D Croitoru
email: mihail.croitoru@u-bordeaux1.fr

A Vagov

A A Shanenko

V M Axt

The Cooper problem in nanoscale: enhancement of the coupling due to confinement

that, no matter how weak the attraction is, two electrons in three-dimensional (3D) space just above the Fermi sea could be bound. In this work we investigate the influence of confinement on the binding energy of a Cooper pair. We show that confinement-induced modification of the Fermi sea results in a significant increase of the binding energy, when the bottom of an energy subband is very close to the Fermi surface.

(Some figures may appear in colour only in the online journal)

As was indicated by Cooper the ground state of the 3D electron gas, i.e. the filled Fermi sphere, is unstable if there is even an infinitesimally small electron-electron attraction near the Fermi surface [1]. Instead, a qualitatively different ground state appears in this system. By considering two electrons above the Fermi sea he found that these electrons are bound, and that the largest binding energy is

E 3D 0 = -2 hω D exp - 2 |g|ν F (1)
for two electrons having opposite momenta (k, -k) and anti-parallel spins (↑, ↓). Here ω D is the Debye frequency, g is a negative constant characterizing the attraction, and ν F is the electron density of states per spin at the chemical potential µ. If electrons having nearly opposite momenta (k, -k + q) are paired (q is the total momentum of a Cooper pair), the binding energy is less than E 0 . The degeneracy of the electron gas (Pauli exclusion principle) is the essential property of the system. Together with the attractive interaction between the particles it is the origin of the Fermi sea instability. It is well known that in the bulk (3D) case a weak attractive potential alone cannot create a bound state [START_REF] Landau | Quantum Mechanics Non-Relativistic Theory 3rd edn[END_REF]. To form a bound state the attraction should be larger than a certain threshold. In contrast, in one-and two-dimensional cases a bound state can be created by an arbitrarily weak attractive potential. It is the filled Fermi sphere at low temperatures that transforms the 3D system to an effectively 2D one. Indeed, the Pauli exclusion principle prohibits the electrons from being scattered inside the Fermi sphere. Besides, the energy transfer between electrons due to the phonon exchange is of the order of hω D , which is much smaller than the Fermi energy, µ. Therefore electrons are effectively confined in a thin shell around the Fermi surface, if hω D µ [START_REF] Zagoskin | Quantum Theory of Many-Body Systems: Techniques and Applications[END_REF]. This result of Cooper made it possible to conclude that the instability of the normal phase, driven by the Pauli exclusion principle and the effective electron-electron correlations is responsible for the occurrence of the superconducting phase [START_REF] Bardeen | [END_REF][START_REF] Schrieffer | Theory of Superconductivity[END_REF].

In recent years, due to remarkable experimental progress [START_REF] Özer | [END_REF][7][8][9][10][11][12][13][14][15][16][17][18], the study of the superconducting properties in nanoscale has attracted a lot of interest [19][20][21][22][23][24][25][26][27][28][29][30][31][32][33]. A characteristic feature of such a system is quantum confinement of charge carriers. It results in important qualitative changes in the characteristics of superconductors with nanoscale dimensions, namely, any superconducting quantity (the superconducting order parameter, the energy gap, the critical temperature [23,24], critical magnetic field [25], critical current [26], coherence length [27]) exhibits quantum-size oscillations with enhancement [9] and reduction [18] trends. In quasi-low-dimensional structures, the spatial distribution of the pair condensate exhibits distinct nonuniform patterns leading to the formation of the Andreev-type states. This significantly modifies the ratio of the energy gap to the critical temperature [32]. The superconducting-to-normal transition driven by a magnetic field becomes a cascade of first-order transitions contrary to the second-order transition found in the mesoscopic regime [25].

The goal of the present paper is to investigate the influence of quantum confinement on the binding energy in the single bound-pair model. The main feature of low-dimensional structures that makes them different from a bulk is the formation of a series of single-electron energy subbands (discrete electron levels in the case of quasi-0D). While the specimen thickness increases (decreases), these subbands move down (up) in energy.

Following Cooper, we consider two quasi-particles above the filled Fermi sphere. We assume translational invariance along one, e.g. the z-direction. We can write the two-particle wavefunction as

(r 1 , r 2 ) = k i,j ψ kij e ik(z 1 -z 2 ) ϕ i (ρ 1 )ϕ i (ρ 2 ), (2)
with i, j the quantum numbers associated with ρ 1 and ρ 2 , respectively, and k the wavevector of the quasi-free electron motion along the z-direction. Here we shall concentrate on the situation where the pair is at rest (q = 0) [START_REF] Zagoskin | Quantum Theory of Many-Body Systems: Techniques and Applications[END_REF]. In this case the scattering phase volume is largest. The two-particle Schrödinger equation reads

(E -Ĥ0 ) (r 1 , r 2 ) = V (r 1 -r 2 ) (r 1 , r 2 ), (3)
where Ĥ0 = -h 2 2m 2 1 + 2 2 . We make an approximation, that the potential V (r 1 -r 2 ) couples only electrons from time-reversed states. We can make this assumption because according to Anderson [34] the pairing of the time-reversed states is a good approximation in bulk in the absence of a magnetic field. For nanoscale systems, the Anderson's approximation is good enough provided that the pairing of two electrons from different single-electron subbands plays a minor role. As shown in [26] this is true when the subband-dependent pairing is much less than the subband energy spacing. We remark that this approximation becomes exact only in the situation of the spatially independent pairing. Our previous numerical investigations [26] have shown that the use of Anderson's recipe is well justified in nanoscale, i.e. corrections were found to be even less than a few per cent. Even for a broken time-reversal symmetry it was shown that Anderson's approximation still provides very good results for a cylindrical superconducting nanowire in a parallel magnetic field [25], because the single-electron wavefunctions remain unaltered in the linear order in the vector potential [25]. Within the Anderson approximation the form of the wavefunction, shown in equation ( 2) can be simplified, namely

(r 1 , r 2 ) = k i ψ ki e ik(z 1 -z 2 ) ϕ i (ρ 1 ) ϕ i (ρ 2 ). (4)
Substituting equation ( 4) into (3) and making projections on both sides of the equation, we obtain the following matrix equation

E -2ε i - h2 k 2 m ψ ki = q j V ki,qj ψ qj , (5)
where V ki,qj = k, i; -k, i|V (r 1 -r 2 ) |q, j; -q, j . This equation is solved with V (r 1 -r 2 ) factorized as

k, i; -k, i|V (r 1 -r 2 ) |q, j; -q, j = - g V ω * ki ω qj , (6)
where the interaction is confined to the area very close to the Fermi surface

|ω ki | 2 = (ξ k,i ), (7)
with

ξ k,i ≡ ε i + h2 k 2
2m -µ, the energy spectrum of an electron in the Fermi sea, measured from the Fermi surface, µ, g is the coupling and V is the volume of the system. Here (ξ ) = 1 if hω D > ξ > 0 and (ξ ) = 0 in other cases. Summing over the quantum numbers k, j, we find the relation

k i ω ki ψ ki = - g V k i ω * ki ω ki E -2ε i - h2 k 2 m × q j ω qj ψ qj . (8)
This is the eigenvalue equation, we have to solve. By measuring the energy of the pair relative to the chemical potential, E = 2µ + E 0 , we rewrite this equation as follows

V g = k i ω * ki ω ki 2ε i + h2 k 2 m -2µ -E 0 . ( 9
)
When the bottom of a parabolic subband passes through the Fermi surface, the density of states ν F increases abruptly, since the density of states per unit length in 1D for one spin projection at the Fermi surface reads as

ν 1D (E) dE = 1 π m 2 h2 1/2 i θ (E -ε i ) (E -ε i ) 1/2 dE, (10)
where θ (Eε i ) is the Heaviside function. Note that the position of the bottom of any subband scales as 1/D 2 , with D the specimen thickness. As an example, in figures 1(a) and (b), the single-electron subbands are schematically plotted versus the wavevector of the quasi-free electron motion along a nanowire. The narrow range of energies for electrons involved in the interaction is highlighted with the patterned boxes. Throughout this paper the energy range [µ, µ + hω D ] is referred to as the Debye window. In figure 1(a) an off-resonant situation is drawn, since the bottom edges of all subbands 1-3 are situated outside of the Debye window. When increasing the thickness of the sample, the subbands move down, and the edge of the third subband now is found in the vicinity of the Fermi level. The resonant regime develops, when a subband enters the Debye window and approaches the chemical potential µ (as shown in figure 1(b)). Size evolution of the system leads to a sequence of peaks in the average density of We can split the sum in equation ( 9) into two terms, namely

V g = k i ξ k,i 2ξ k,i -E 0 + k ξ k,R 2ξ k,R -E 0 , (11)
where the prime means the absence of resonant subbands in the summation. Here it is supposed for simplicity that only one resonant subband is present in the spectrum that is labeled by i = R. Considering i as a continuous variable and making an approximation that the contribution from the non-resonant subbands to the density of states weakly depends on the energy in the Debye window, i.e. ν(E) = ν(µ), since for the non-resonant subbands hω D µ i ≡ µε i , we can estimate the first term in equation (11) as

F b = g hω D 0 ν(E) dξ 2ξ -E 0 = g ν F 2 ln E 0 -2 hω D E 0 , (12)
where ν F is the density of states per unit volume at energy µ.

In the evaluation of the second term in equation ( 11), which we denote as F c , we make use of equation (10) for the density of states per unit length in 1D. Then F c acquires the following form

F c = g πL 2 m 2 h2 1 2 hω D 0 dξ 2ξ -E 0 1 (ξ -ξ R ) 1/2 , ( 13
)
where ξ R = ε R -µ = -µ R is the energy of the bottom of the resonant subband with respect to the chemical potential and L is the nanowire thickness. ξ R is negative because we consider the regime where the bottom of the resonant subband is below the chemical potential. Performing the integration we get

F c = g πL 2 m h2 ζ 1 2 × arctanh - 2ξ R ζ 1 2 -arctanh 2 hω D -2ξ R ζ 1 2 (14) with ζ = E 0 -2ξ R = E 0 + 2µ R .
If the bottom edge of the subband R is situated very close to the Fermi level, we can further approximate ξ R = 0 and

F c = g πL 2 m h2 E 0 1 2 arctanh 2 hω D E 0 1 2 . ( 15
)
Since usually hω D |E 0 |, and making use of lim

hω D |E 0 | arctanh - 2 hω D E 0 1 2 = iπ 2 , (16)
we get the total condition in the form

1 g = ν F 2 ln E 0 -2 hω D E 0 + 1 L 2 - m 4 h2 E 0 1 2 . (17)
Of course if there are M resonance subbands, with bottom edges in the vicinity of the Fermi surface the second term should be multiplied by M.

If the system is in the bulk limit, L → ∞, then the contribution of the second term is negligible and we obtain the Cooper expression for the binding energy (1). The reduction of dimensions of the system when approaching the nanoscale regime splits the conduction band into a multiple of subbands. Since confinement is weak, the quantization energy is small and the number of subbands with edges in the Debye window is large (see figure 2(c)). There is always a finite and large probability of finding the edge of a subband very close to the Fermi surface. This subband we call resonant. Since the number of subbands is large, the contribution of the resonant subband to the total density of states is small and can be considered as a perturbation. For the zero-order perturbation F c = 0 and this leads to the bulk result of equation (1). For the first-order perturbation we obtain the Cooper result modified due to the influence of the resonant subband

E 1D→3D 0 = -2 hω D exp - 2 | g|ν F , ( 18
)
where

g -1 = g -1 - 1 hν F L 2 m 2 hω D 1 2 exp 1 |g|ν F . ( 19
)
Hence the binding energy of the pair increases exponentially as L decreases. This is the origin of the sharp resonances.

In the case of strong confinement, the resonant subbands' contribution to the density of states is significant, even determinative (see figure 2(a)). When the edge of the resonant subband touches the Fermi surface the main contribution to the density of states comes from the second term in equation (17). The contribution from the other subbands with bottom edges below the Fermi surface can be viewed as a perturbation. In the zero-order approximation equation ( 17) reads as

1 g = 1 L 2 - m 4 h2 E 0 1 2 , (20)
which leads to the well-known 1D result for the energy of the bound state in the attractive δ-function potential with the coupling constant g/L 2 (we have the extra factor 1/2 due to the reduced mass m/2)

E 1D 0 = - mg 2 4 h2 L 4 . (21)
Making the same iterative procedure as above, we get the Cooper pair binding energy in the first-order perturbation

E 3D→1D 0 = - mg 2 4 h2 4 1 -|g|ν F 2 ln 8h 3 ω D L 4 mg 2 2 . (22)
The model presented here is oversimplified. For example it cannot describe a situation of moderate confinement, as presented in figure 2(b). Nevertheless it enables us to investigate the influence of the system's dimensionality on the binding energy of the Cooper pair in the presence of the Fermi sea. Notice that the binding energy obtained in the zero-order approximation, equation (21), is just the exact binding energy of two particles in a vacuum for the case of a 1D interaction model with V(x) = -V 0 δ(x). Therefore, quantum confinement in the resonant-subband regime opens a channel for the real (molecule) binding. In other words, the presence of the resonant subband very near to the Fermi surface reconfigures the scattering phase space and forms an effective low-dimensional situation. The physics behind can be understood as follows. Two additional electrons can scatter via their interaction to all states of the same total momentum that are unoccupied by other electrons. The correlation between these electrons depends on the number of states in the Debye window above the Fermi sea available to form a pair [29]. When a subband approaches the Fermi level the number of states available for scattering strongly increases. This results in an increase of the binding energy. Since the Fermi sea in our model is frozen, i.e. there is no correlation between the electrons in the Fermi sea and two additional electrons, these two electrons can scatter only from the present state to all other empty states, but there is no in-scattering flow to these states. Besides, there is no possibility for these two electrons to scatter below the Fermi surface.

So far we have studied pairing in the situation of ξ R < 0. Let us consider another regime when the chemical potential is situated slightly below the bottom edge of the empty subband. We can also call this subband resonant and split the sum in equation ( 9), resulting in equation (11), because the contribution to the averaged density of states from this subband is larger than from the others. Hence, equation (11) reads as

1 = g 2πL 2 2m h2 1 2 hω D ξ R dξ 2ξ -E 0 1 (ξ -ξ R ) 1/2 + F b . (23)
In the zero-order perturbation, the integration of this equation gives

1 = g 2πL 2 2m h2 1 2 arctanh hω D -ξ R (ξ R -E 0 /2) 1 2 (ξ R -E 0 /2) 1/2 . ( 24
)
Assuming that hω D |E| and hω D ξ R we find that

E 0 = 2ξ R - mg 2 4 h2 L 4 . (25)
The first term here is positive and the second is negative. While the bottom edge of the resonant subband moves up, the binding energy in this approximation decreases and finally it becomes zero, when 2ξ R = mg 2 4h 2 L 4 , meaning that the contribution of the resonant subband is no longer dominant and we have to take into account the contributions from the other subbands. In other words, we are back to the bulk solution for the Cooper pair problem: the resonance in the binding energy disappears. We remark that the same occurs when the bottom of the resonant subband goes far below the Fermi surface.

It is known that the spatial extent of the pair in the 3D Cooper problem is ξ 3D ≈ 1.15 hv F /E 0 , which is ∼10 4 Å. In the resonant regime, the pair binding is enhanced and it is expected that this enhancement should modify the size of the Cooper pair. Therefore, it is of interest to investigate the influence of quantum confinement on the size of a Cooper pair. Obviously, the size of the pair in the confined directions is determined by the size of the confinement. The size of the Cooper pair along the non-confined direction may be defined as

ξ 2 C = | (ρ 1 , ρ 2 , z)| 2 z 2 dz | (ρ 1 , ρ 2 , z)| 2 dz = i ϕ i (ρ 1 )ϕ i (ρ 2 ) ∂ ∂k ψ ki 2 dk i ϕ i (ρ 1 )ϕ i (ρ 2 ) |ψ ki | 2 dk. ( 26
)
If we neglect the contribution from the non-resonant subbands and consider in the summation only the resonant one, substituting equation ( 4) here we obtain

ξ 2 C = hω D 0 dξ (2ξ -E 0 ) 4 (ξ -ξ R ) 1/2 hω D 0 dξ (2ξ -E 0 ) 2 1 (ξ -ξ R ) 1/2 , ( 27
)
where we have taken into account that

(E 0 -2ξ R )ψ ki = q j
V ki,qj ψ qj = -λω * k,i q j ω q,j = const, (28) and performing the integration under the condition hω D → ∞ we get

ξ 2 C = h2 8mE 0 . ( 29
)
Since in the resonance subband the binding energy E 0 is large, the size of the resonant Cooper pair becomes small. We can say that the Cooper pair becomes like a point boson dimer. In summary, we solved the Cooper problem of two fermions attracting each other pairwise via a 3D potential in the presence of quantum confinement. We found that in the strong confinement regime the behavior of the binding energy is qualitatively different from that in the bulk. In the presence of quantum confinement the band of single-electron states in clean nanosample is split up in a series of subbands. While the specimen thickness increases (decreases), these subbands move down (up) in energy. Each time when the bottom of a parabolic subband passes through the Fermi surface, the binding energy of a resonant Cooper pair strongly increases. The resonant Cooper pair becomes like a point dimer.

Figure 1 .

 1 Figure 1. (a) Scheme of the single-electron subbands versus the wavevector of the motion parallel to the specimen for the off-resonant case. (b) The same as in panel (a) but for the resonant thickness. (c) The 1D density of states at the Fermi level versus the wire thickness D.

Figure 2 .

 2 Figure 2. Scheme of the single-electron subbands versus the wavevector of the motion parallel to the specimen for (a) an ultra-small sample, (b) a sample with moderate confinement, (c) a weakly confined sample.

Acknowledgment MDC acknowledges support by the European Community under the Marie Curie IEF Action (Grant Agreement No. PIEF-GA-2009-235486-ScQSR).


David Poirier-Quinot
email: david.poirier-quinot@sorbonne-universite.fr

Brian F G Katz

David Poirier

Impact of HRTF individualization on player performance in a VR shooter game II

. Results suggest that the impact of the HRTF on participant performance (speed and movement efficiency) depends both on participant sensitivity and HRTF presentation order.

Introduction

"Binaural hearing" refers to the capability of integrating information from the two ears to perceive a sound in three-dimensional space (azimuth, elevation, and distance). Psychophysical studies have shown that various mechanisms are involved in the human auditory system for sound localization [START_REF] Blauert | Spatial Hearing: The Psychophysics of Human Sound Localization[END_REF]. To infer the angular direction of a sound source, these mechanisms rely on directiondependent audio cues, resulting from the propagation of an acoustic wave from the source to both ears. Using digital signal processing, these cues can be applied to any audio input to simulate a sound object at a virtual position in a listener's 3D auditory space (experienced over headphones). The set of these direction-dependent cues for a given person is typically referred to as a Head Related Transfer Function (HRTF).

HRTFs are individual, directly resulting from the interactions between a person's morphology and an impinging acoustic wave during its propagation around the head [START_REF] Blauert | Spatial Hearing: The Psychophysics of Human Sound Localization[END_REF]. Individuals listening to a binaural audio scene rendered using their own HRTF will perceive each of its components with more spatial precision than those presented with a random HRTF set [START_REF] Begault | Direct comparison of the impact of head tracking, reverberation, and individualized headrelated transfer functions on the spatial perception of a virtual speech source[END_REF]. Various methods have been proposed to select a "best-match" HRTF from an existing database [START_REF] Seeber | Subjective selection of non-individual head-related transfer functions[END_REF][START_REF] Zotkin | HRTF personalization using anthropometric measurements[END_REF], as measuring an in-dividual's actual HRTF is a demanding operation [START_REF] Carpentier | Measurement of a head-related transfer function database with high spatial resolution[END_REF]. This process of selection and use of a best-match HRTF for binaural synthesis is here referred to as HRTF individualization.

This study is part of a larger research project aiming to characterize the impact of binaural rendering quality in the context of different Virtual Reality (VR) application contexts. This study focuses on the impact of HRTF individualization on performance in the context of a VR shooter game. While the core of the gameplay is built around an audio-visual localization task, it extends the existing literature [START_REF] Wenzel | Localization using nonindividualized head-related transfer functions[END_REF][START_REF] Xu | Individualization of head-related transfer function for threedimensional virtual auditory display: a review[END_REF][START_REF] Härmä | Personalization of headphone spatialization based on the relative localization error in an auditory gaming interface[END_REF][START_REF] Mehra | Comparison of localization performance with individualized and non-individualized headrelated transfer functions for dynamic listeners[END_REF] in that the final experience is truly a game, where participants are placed under increasingly difficult time and performance constraints.

The underlying hypothesis of this study is that using individualized HRTFs will result in an increase in participant performance (reaction time, movement efficiency), more so as the overall game dynamic (enemy spawn interval, flight speed, etc.) increases.

Experimental design

The experiment consisted of two sequential parts. The objective of Part 1 was to identify the best and worst match HRTFs from a subset of 7 for each participant. Part 2 was the VR shooter game. A total of 20 participants undertook the experiment (8 female, mean age 31.8 ± 11.5 years).

HRTF classification

The HRTF subset database was assembled from the LISTEN database, defined from a "perceptually orthogonal" optimized HRTF collection [START_REF] Katz | Perceptually based head-related transfer function database optimization[END_REF]. Per-participant best and worst match HRTF sets were selected based on the method elaborated in [START_REF] Andreopoulou | Subjective HRTF evaluations for obtaining global similarity metrics of assessors and assessees[END_REF], establishing a classification based on perceptual-space distance between a spatialized audio trajectory and a described reference. Two trajectories were presented: horizontal plane (12 angles [0°:30°:330°] ) and median plane (19 angles [-45°:15°:225°]), as illustrated in Figure 1.

Each audio trajectory was generated with the 7 HRTFs from the subset. Participants were instructed to rate the 7 resulting versions of each trajectory on a fixed 9-point scale. They were encouraged to distribute their notations on that scale, and required to indicate at least one best [START_REF] Xu | Individualization of head-related transfer function for threedimensional virtual auditory display: a review[END_REF] and one worst (1) match. Both median and horizontal rating sessions were repeated 3 times, to take into account HRTF rating consistency [START_REF] Andreopoulou | Investigation on Subjective HRTF Rating Repeatability[END_REF].

Participants completed the study in a listening booth, ambient noise level < 30 dBA using Sennheiser HD600 headphones and RME Fireface UC audio interface.

VR Shooter Game

During the VR shooter game, participants were equipped with an Oculus CV1 Head Mounted Display (HMD), a pair of headphones (those of the CV1), and a pair of hand tracked devices (Oculus Touch controllers).

The game started with participants immersed in a virtual scene, standing on a 0.5 m radius platform mounted on a pole at the center of a 20 m radius spherical structure. Enemy targets could "spawn" from any of the 29 evenly distributed holes in the structure, flying in straight lines towards the participant until collision, either with a bullet or the participant. Participants were instructed to shoot at the incoming targets using a pair of hand-held blasters, the avatar representations of the hand tracked devices in the virtual scene, destroying as many as possible, as fast as possible, in the given time limit.

Enemy targets emitted specific event-based sounds for: spawning, launching, flight, and collision. All sounds were spatialized using the Anaglyph binaural audio engine v0.9 [START_REF] Poirier-Quinot | The Anaglyph binaural audio engine[END_REF]. Anechoic conditions were employed, no room effect was included to keep the study's focus on HRTF effects.

The game was designed using the Unity v2017.3.0 game engine with modelled assets designed in Blender v2.79. The OSC (Open Sound Control) protocol [START_REF] Wright | Open Sound Control: an enabling technology for musical networking[END_REF] was used for communications between Unity and the Anaglyph engine running as a VST in Cycling'74 Max v7.3. Figure 2 depicts the game setup and the VR scene.

A video extract showing a game session is available 1,2 .

A short training session introduced the controls and the difficulty level mechanism, implemented so that the overall game dynamic (enemy spawn interval, flight speed, etc.) increased as the game progressed and participant's skill (in-game "level" in Table 1) improved. The 5 min game was played in six sessions, alternating between participant's best and worst match HRTF. Best/worst match presentation order was evenly balanced among participants, resulting in two groups. To avoid fatigue on this rather demanding task (see video extract), participants played sessions S1-2, directly following the HRTF classification task, followed by a week pause interval, then sessions S3-6.

Results

Preliminary elements of the results were previously presented [1]. That previous study was based on the current experimental design and protocol, involving 30

1 https://youtu.be/q6muds1qW-w 2 https://www.youtube.com/c/ LAMSorbonneUniversité participants for just the first two game sessions. The present study extends these results to a six session game so as to investigate the impact of HRTF individualization on the "long" term, to dissociate it from the game learning effect that may have affected the results of the first part of this study. All of the 20 participants involved in this study took part in the first experiment, each undertaking four more sessions approximately a week after they took the first two.

Performance assessment was based on three metrics, calculated from sessions logs: in-game level, spawnspot reaction time and spawn-spot travelled angular distance. The in-game level was related to the number of enemies destroyed versus those that hit participants: increasing one unit for every three consecutive kills, decreasing one unit for every two consecutive fails. Participant spawn-spot reaction time corresponded to the time interval between the spawn of a target and its entering the visual field of view, defined as a 50°cone centered around the current forward view axis. The event of seeing the target, rather than destroying it, was chosen so as to remove the impact of skill at aiming and destroying targets from the analysis (the task of targeting being independent of the acoustic rendering quality). Targets visible upon spawn (spawned in the current field of view) were discarded from spawn-spot reaction time analysis. For targets that never entered participants field of view, colliding with or shot by the participants, the spawn-to-collision time was used as the spawn-spot reaction time. Such targets, shot without ever being visible, represented less than 3% of the total number of target spawned.

The associated spawn-spot travelled angular distance metric corresponds to the angular distance traversed by the participant's head from target spawn to target spot sub-events. This last metric represented movement efficiency, and served to differentiate between participants using binaural cues to localize targets and those randomly looking around [START_REF] Katz | The use of 3D-audio in a multi-modal teleoperation platform for remote driving/supervision[END_REF].

Table 1 summarizes the independent and dependent variables of the experimental protocol. Result significance was assessed using a Wilcoxon signed rank test (p-value threshold of 0.05) as all compared pairedsample distributions proved to follow a non-normal (skewed) distribution.

HRTF Classification Results

Results of the HRTF ratings of Part 1 are summarized in Figure 3, focusing on the scores obtained by each participant's best and worst HRTF match for both trajectories.

Participants were consistent in their classification with regards to these extrema. As audio sources in Part 2 of the experiment were to arrive from all directions, an average best-and worst-HRTF match across trajectories was established for each participant. The rating statistics for selected best-and worst-HRTF matches are shown in Fig. 3c. Most participants proved consistent in their ratings to clearly distinguish between best and worst match for both horizontal and median trajectories. As also shown in [START_REF] Andreopoulou | Subjective HRTF evaluations for obtaining global similarity metrics of assessors and assessees[END_REF], participant HRTF ratings for the horizontal and the median trajectories were not correlated. This explains the observed decrease in rating values for the trajectory mean results (Fig. 3c) as compared to the individual horizontal and median plane trajectory ratings (resp. Fig. 3a and Fig. 3b). For almost all participants, average-best and average-worst HRTF scores remained sufficiently distinct to distinguish both populations.

VR Shooter Game Results

Result analysis was subdivided into session wise and event wise analysis. Session wise analysis concerns participant mean results across sessions. Event wise analysis decomposes each session in events, an event being defined as an enemy target {spawn, launch, flight, collision} sequence. Events related to targets spawned in participant's field of view (see Section 3) have been discarded from analysis. Events where the target was not destroyed and ended up colliding with the participant are included in the analysis, adding both total event time and angular distance to participant's results.

Out of 20 participants, 10 started S1 with their worstmatch HRTF (group 1), 10 with their best-match HRTF (group 2). The analysis below concerns 20 × 6 = 120 sessions for a total of 14430 events (average of 120.3 ± 18.3 events per session, since game dynamics vary with participant in-game level). plateau from S4 on (from 14.6 for S1 to 18.1 for S4-6). Mean spawn-spot angular distance traversed and response time (both averaged per participant session) likewise significantly diminished from S1 to S4 (from 147.7°and 1.37 s for S1, to 125.0°and 1.27 s for S4-5-6 resp.). These results highlight the game task learning effect independent of the HRTF condition.

Statistical analysis across participants

Repeating the same analysis as a function of HRTF across both participant groups (Figure 4b) showed no significant impact of the HRTF on these metrics. These results would suggest at first glance that there was no benefit to using individual HRTFs in the VR shooter game.

Statistical analysis by event

Combining the results of each event for all participants, the event-wise analysis also reflects the significant impact of game learning on performance, for spawnspot angular distance traversed (from 147.2°for S1 to 125.1°for S6) and reaction time (from 1.35 s for S1 to 1.26 s for S6). As for mean session results in Section 3.2.1, event-wise statistics show no significant impact of the HRTF quality on either reaction time or angular distance traversed. This result differs from those reported in [1], limited to the first two sessions of the game. HRTF quality still has a significant impact on both metrics between the first two sessions for these 20 participants, yet none for all sections combined, nor for any S3 to S6 combination.

HRTF×Session interaction analysis

Subsequent analysis examines the interaction effect of HRTF presentation order across sessions. The eventwise mean and 95% confidence intervals for reaction time and angular distance traversed are shown in Figure 5. As seen in [1, Figure 5], a clear interaction effect can be observed on event-wise reaction time between S1 and S2, only for the group of participants going from worst to best HRTF. This suggested at the time that "a worst match HRTF negated the benefit that should result from training". The extended results of sessions S3-S6 show no further impact of HRTF quality, regardless of the presentation order.

Participants undertook S1-2 immediately after the HRTF classification task, S3-6 took place at least a week after that. The difference in performance improvement from S1 to S2 between both groups not being repeated in S3-6 could suggest either a difference in game related learning capability between the two groups, or a desensitisation to HRTF quality as acquired after the classification task. The week gap between S2 and S3 did not result in any clear performance drop, suggesting that participants did not "unlearn" the game between both sessions.

As seen in [1], no significant difference was observed between the performance of both groups in S1. In the long run, participants from the first group (started with worst) significantly out-performed those from the second group (average results from S3 to S6: 127.3°and 1.24 s versus 133.0°and 1.32 s), overall showing a steeper learning curve regardless of the HRTF quality. This difference in learning ability between both groups weakens the hypothesis made in [1] on the worst match HRTF negating the benefit that should result from training, which indicates rather an overall performance bias for group 1 participants with respect to group 2.

HRTF×Spawn elevation interaction analysis

Figure 6 presents event-wise angular distance traversed for target spawn elevation and HRTF. Spawn elevation regions R2-4 resulted in significantly lower distributions than R1 and R5-6 (121.7°versus 156.4°) due to participants resting head orientation (facing R3 or R2 most of the time), typical search pattern (azimuth search followed by elevation search), and up/down confusions induced by binaural rendering. As could be expected, participants had difficulties locating targets when spawned from the R6 region, significantly more so when using their worst HRTF (mean angular distance of 158.1°versus 176.0°for best and worst resp.). On average, participants only managed to locate and destroy half of the targets spawned from R6 (54% out of 281 spawns with best versus 44% out of 282 spawns with worst). No significant impact of HRTF quality on angular distance traversed was observed for any other spawn elevation region. No significant region×HRTF interaction was observed for participant reaction time.

Post hoc analysis of participant sensitivity to individualized HRTF

Correlation analysis of per-participant mean angular dist. versus reaction time shows poor results as a function of best and worst match HRTF sessions (r = 0.30 Radius of concentric white circles are proportional to mean variance of best and worst match HRTF ratings (see Fig. 3). Numbers are participants game performance rank, based on the best average level across all sessions. Blue and red patches are ellipse fits around blue and red clusters resp. using a least squares criterion.

and r = 0.26 for best and worst resp., for combined spawn elevation regions R1 and R5-6). In contrast, examination of the per-participant mean differences between best and worst match HRTF angular dist. versus reaction time, again for spawn elevation regions R1 and R5-6, shows a strong correlation of r = 0.89, see Fig. 7.

Clustering of participant performance based on HRTF presentation order is highlighted by enveloping ellipses. This clustering indicates that most participants who performed best with their best match HRTF (resp. best with worst match) started S1 their worst match (resp. with their best match). This result suggests that HRTF presentation order had an impact on HRTF-wise performance, regardless of participant specific HRTF ratings. The HRTF presented in session S2 seems more likely to be "efficient" than the one presented during the first session, where learning the game takes precedence over focus on audio localization. This result, and its interpretation, would deserve a more thorough study, as it could simply be related to a nonuniform distribution of participants ability among both groups as observed in Section 3.2.2. Worth noticing, 4 of the top 5 participants (mean level based ranking) started with their worst HRTF and are located in the upper-right quadrant, i.e. performed best with their best HRTF.

The clustering of Figure 7 strongly suggests that some participants are more sensitive to HRTF subjective quality than others, or that at least the classification resulting from the HRTF rating was more appropriate for these participants for the task at hand. Out of 20 participants, HRTF quality proved to have a significant impact on event-wise angular distance traversed for 4 of them, 3 who performed best with their best match HRTF, 1 with the worst match ({1, 4, 14} vs{6} in Fig. 7). For none of the participants did HRTF quality show a significant impact on spawn-spot reaction time.

No clear correlation could be observed between point coordinates and radius, i.e. between best-worst HRTF performance difference and best-worst HRTF ratings differences. Likewise for HRTF rating variance (points inner circle size, Fig. 7). This result does not support the "creation of an overall metric to rate participants affinity with binaural hearing [based on participant results to the HRTF classification task]" as suggested in [1].

Conclusion

We have presented the results of an experiment designed to assess the impact of individualized binaural rendering on player performance in the context of a VR "shooter game". Participants performed six game sessions, alternatively using their best-and worst-match HRTF, extending the two-session game published in [1] to further dissociate the impact of the learning effect and that of the HRTF on participants performance. During the game, participants had to locate and shoot at successive enemy targets approaching from random directions within a sphere.

Results indicate that the use of a best-match HRTF did not improve overall participants performance regarding the time they needed to localize the targets and the angular distance they travelled before doing so. An analysis focused on spawn region however revealed that angular distance traversed significantly decreased when participants used their best match HRTF for the top-most R6 region of Fig. 6 (mean of 158.1°versus 176.0°for best and worst resp.). Targets spawned in this region were also more often spotted before collision by subjects using their best-match HRTF (54% versus 44%).

Participant clustering in Fig. 7 indicated that 3 out of 20 participants were significantly more efficient in term of angular distance traversed to locate targets with their best HRTF. This result strengthen the hypothesis formulated in [1] that "the benefits of HRTF individualization for a sub-group of "aware listeners" would exceed those of the average participants". While these aware listeners accounted for less than a fifth of the tested participants, it may be that another HRTF selection method, or a different game design, could increase this ratio. A potential impact of the HRTF presentation order on participants sensitivity to best and worst match HRTF was inferred from the results of this clustering (see Section 3.2.5). This last observation would require further study before any assertion can be made.

A last extension of this study is currently being conducted, adding a control group for reference to characterize the learning effect for a fixed set of best versus worst match HRTF during the early stages of the game. A final and more thorough statistical analysis on participant results will then be conducted, to conclude on the impact of HRTF quality on participant performance for the task at hand.

Fig. 1 :

 1 Fig. 1: Trajectory descriptions for HRTF quality ratings: horizontal (left) and median (right) plane trajectories indicating the start/stop position and trajectory direction (-→).

Fig. 2 :

 2 Fig. 2: Game scene overview: (upper-left) overall view of the virtual environment, (upper-right) focus on the platform atop which participants stand during the game, (lower-left) participant in the VR room, (lower-middle) virtual environment during gameplay, and (lower-right) in-game screenshot.

Fig. 3 :

 3 Fig. 3: Results of the HRTF classification task for all participants. The reported rating values correspond to the average normalized rank given by participants to their determined best/worst HRTF. A value of 0 (resp. 1) indicates that the HRTF was always rated as the least (resp. most) representative of the described trajectory across the 3 rating repetitions. Scores for the best and worst HRTF matches for (a) horizontal and (b) median trajectories. (c) Combined trajectory mean score results for the determined best and worst match HRTF. Error bars indicate the variance of participant ratings for best and worst HRTF.

Figure 4 Fig. 4 :

 44 Figure 4 illustrates the evolution of participants ingame level across sessions. Mean in-game level significantly increased across sessions, until it reached a

Fig. 5 :

 5 Fig. 5: Agregated (a) event-wise reaction time and (b) normalized angular distance (mean and 95% CI) across sessions and HRTF for both groups of participants.

Fig. 6 :

 6 Fig. 6: angular distance across target absolute spawn region, separating worst best match HRTF. R i represent the spawn region. All spawn positions of a given region share a common elevation.

Fig. 7 :

 7 Fig. 7: Participants clustering based on the difference in event-wise performance using best and worst match HRTF for spawn elevation regions R1 and R5-6 combined. Each point represent the performance of a given participant, abscissa and ordinate values are differences in mean response time and angular dist. between best and worst match HRTF sessions respectively. Points in the top-right section (resp. bottom-left) represent participants who best performed with their best (resp. worst) match HRTF. Blue circles are participants who started the game with worst match HRTF, red for those who started with their best. The radius of each circle is proportional to the difference between best and worst match HRTF scores obtained during rating task.Radius of concentric white circles are proportional to mean variance of best and worst match HRTF ratings (see Fig.3). Numbers are participants game performance rank, based on the best average level across all sessions. Blue and red patches are ellipse fits around blue and red clusters resp. using a least squares criterion.

Table 1 :

 1 Independent and dependent variables of the experimental protocol.

	Independent variables
	Participant ID	20	random variable
	HRTF ID	7	best, worst
	Session ID	6	S1, S2, ..., S6
	Spawn region	6	R1, R2, ..., R6
	Dependent variables
	angular distance	event-wise	raw and norm
	time	event-wise	raw and norm
	mean level	session-wise	raw and norm

Acknowledgements

This work was funded in part through a fundamental research collaboration partnership between Sorbonne Université, CNRS, Institut ∂ 'Alembert and Oculus VR, LLC.