X Zhou

V Schmitt

P Bertet

D Vion

W Wustmann

V Shumeiko

D Esteve

High-gain weakly nonlinear flux-modulated Josephson parametric amplifier using a SQUID-array

We have developed and measured a high-gain quantum-limited microwave parametric amplifier based on a superconducting lumped LC resonator with the inductor L including an array of 8 superconducting quantum interference devices (SQUIDs). This amplifier is parametrically pumped by modulating the flux threading the SQUIDs at twice the resonator frequency. Around 5 GHz, a maximum gain of 31 dB, a product amplitude-gain × bandwidth above 60 MHz, and a 1 dB compression point of -123 dBm at 20 dB gain are obtained in the non-degenerate mode of operation. Phase sensitive amplification-deamplification is also measured in the degenerate mode and yields a maximum gain of 37 dB. The compression point obtained is 18 dB above what would be obtained with a single SQUID of the same inductance, due to the smaller nonlinearity of the SQUID array.

Although superconducting parametric amplifiers based on Josephson junctions are known and understood for decades [START_REF] Barone | Physics and applications of the Josephson effect[END_REF][START_REF] Yurke | [END_REF], they have recently received an increased attention [3] because of their ability to measure single quantum objects and engineer quantum fluctuations of a microwave field. They are extensively used to readout superconducting quantum bits [4,6,7] or mechanical resonators [8] at or near the quantum limit, i.e with minimum back-action imposed by quantum mechanics for the given amount of information taken on the system. They permitted for instance the measurement of quantum trajectories [4] and the implementation of quantum feedback schemes [5,9]. In the field of quantum microwaves, they are also used to squeeze quantum noise and produce itinerant squeezed states for encoding quantum information [10,11] or demonstrating fundamental effects like the reduction of the radiative decay of an artificial atom [12].

Compared to the noisier high electron mobility transistor (HEMT) based amplifiers, these Josephson parametric amplifiers (JPA) suffer from limited bandwidth and from gain saturation at extremely low input power. A strong effort is thus made to increase the bandwidth and to mitigate saturation of JPAs by varying their design and mode of operation [13][14][15][16]. In all cases, parametric amplification of a signal at angular frequency ω S occurs by transfer of energy from a pump at frequency ω P to the signal and to a complementary idler frequency ω I . For amplifiers based on resonators, one distinguishes the case of intrinsically nonlinear resonators with bare frequency ω R that are pumped at ω P ω S ω R directly on their signal line, from the (possibly linear) resonators whose frequency ω R is parametrically modulated with a pump tone at ω P 2ω S 2ω R on a dedicated line separated from the signal port. In the first case the intrinsic nonlinearity of the resonator is usually obtained by implementing all or part of its inductance by Josephson junctions (or superconducting weak links). The pumping at ω P at sufficiently high amplitude modulates this nonlinear inductance at 2ω P , and is responsible for a 4wave mixing such that 2ω P = ω S + ω I . In the present work, we are interested in the second case [17,18], for which the nonlinearity is due to an externally imposed parametric modulation of the frequency and is responsible for a 3-wave mixing such that ω P = ω S + ω I . The interest of this 3 wave mixing is that no pump mode propagates along the input and output signal lines and can blind a detector or spoil a squeezed field, at a close frequency. In practice, the true parametric modulation is usually obtained by embedding in the resonator inductance one or several SQUIDs, the Josephson inductance of which is modulated by an ac magnetic flux. The nonlinearity of the resonator inherited from the SQUID(s) is in this case an unwanted feature, which leads to saturation of amplification, and should thus be kept low. So besides the advantage of getting rid of the pump along the signal lines, a truly parametrically pumped amplifier can also be made more robust against saturation by reducing its nonlinearity without having to pump it more strongly. In this work, we test this idea and demonstrate a weakly nonlinear JPA with high gain, made of a lumped LC resonator with the inductor L terminated by a SQUID array. The manuscript first summarizes the theoretical description of such a JPA, then describes the device implemented and its characterization setup, and finally presents the experimental data and a comparison between measured and calculated gain, bandwidth, and saturation.

I. THEORETICAL SUMMARY AND DESIGN CHOICES

The specificity of the JPA presented here (pure 3 wavemixing with spurious nonlinearity) makes the standard classical description of parametric amplifiers [START_REF] Yurke | [END_REF] not directly applicable to it. This is why a comprehensive theoretical summary is given here, based on the theoretical work [21] (note that a similar theoretical treatment can be found in [19]). The equivalent circuit of the JPA is shown in the bottom-right corner of Fig. 1b. For a DC flux bias Φ DC , a parametric modulation δL A cos(ω P t) of its array inductance L A , and a microwave input signal V S cos(ω S t + χ), the JPA equation of motion at the lowest nonlinear order in phase ϕ = 1/ϕ 0 ´V dt across the total inductance

L = L g + L A (Φ DC ) or the capacitance C R is ϕ +2Γ a ϕ +ω 2 R [1+a P (Φ AC ) cos(ω P t)]ϕ+α 1 ϕ 3 = ϕ S cos(ω S t+χ), (1)
where

ϕ 0 = /2e = Φ 0 /2π is the reduced flux quantum, ω R = 2πf R = 1/ L(C R + C c
) the frequency of the resonator, ω P = 2πf P the pumping frequency, Γ a = ω R /2Q its amplitude decay rate, a P = δL A (Φ AC )/L the relative pumping amplitude, α 1 = -ω 2 R p 3 /2N 2 the Josephson nonlinearity coefficient with N the number of SQUIDs and p = L J /L the so-called participation ratio of the total Josephson inductance L J to the total inductance L; ϕ S is the drive amplitude proportional to V S . Taking into account the finite ratio β of each SQUID loop inductance to the inductance L J1 of a single junction, the SQUID array inductance is [20] and x = πΦ/Φ 0 .

L A = N L J1 β/4 + L J with L J = N L J1 / 2 cos (x) + β/2 sin 2 (x)
Equation (1) contains the parametric nonlinearity cos(ω P t)ϕ and the intrinsic Josephson nonlinearity α 1 ϕ 3 mentioned in the introduction. We rewrite it in the frame rotating at ω P /2 using the slow complex internal amplitude A(t) defined by ϕ(t) = A(t) √ Z R e -iωPt/2 + cc, as well as the constant complex amplitude B S of the input signal and the slow output amplitude C(t) related to the input and output voltages V in,out by V in (t) = iϕ 0 √ Z 0 ω S B S e -iωSt /2 + cc and V out (t) = iϕ 0 √ Z 0 ω S C(t)e -iωSt /2 + cc. Here, Z 0 is the impedance of the line, Z R = L/(C R + C c ) the characteristic impedance of the resonator, and cc stands for the complex conjugate of the previous term. Neglecting fast oscillating terms (rotating wave approximation), one obtains [21] iA

+ (iΓ a + δ + α A 2 )A + A * = √ 2Γ a B S e -i∆t C = -i √ 2Γ a A + B S , (2)
where δ = ω P /2 -ω R is the pump to resonator detuning, ∆ = ω S -ω P /2 the signal to pump detuning, α = -p 3 Z R ω R /16N 2 the new nonlinear coefficient, and = ω R a P = 2ω R κΦ AC /Φ 0 the pumping strength with κ ∝ pQ the relative frequency change per flux quantum deduced from the slope of the modulation curve ω R (Φ DC ).

Although the most general stationary solution of Eq. ( 2) is a sum k∈Z A k e -ik∆t of all harmonics at frequencies ω P /2 + k∆, only the signal A S = A 1 and the idler A I = A -1 contributions happen to be non negligible at not too high pumping strength. In this case, they obey

           δ r + α r |A S | 2 + 2 |A I | 2 + ∆ r + i A S + r A * I = 2/Γ a B S δ r + α r 2 |A S | 2 + |A I | 2 -∆ r + i A I + r A * S = 0 C S = -i √ 2Γ a A S + B S C I = -i √ 2Γ a A I (3) with δ r = δ/Γ a , ∆ r = ∆/Γ a , α r = α/Γ a = -p 3 Z R Q/8N
2 and r = /Γ a the dimensionless detunings, nonlinear coefficient, and pumping strength, respectively.

Our goal is to make the nonlinearity α r as small as possible and benefit from the linear signal and idler complex gains given by system (3) when α r = 0, i.e.

G S = C S B S = δ 2 r -∆ 2 r -1-2 r -2iδr δ 2 r -∆ 2 r +1-2 r -2i∆r G I = C I B S = 2i r δ 2 r -∆ 2 r +1-2 r +2i∆r e i2χ , (4)
yielding the power gains

|G S | 2 = 1 + |G I | 2 = 1 + 4 2 r [1 -2 r + δ 2 r -∆ 2 r ] 2 + 4∆ 2 r . (5)
In the degenerate case corresponding to ∆ = 0, the signal power gain becomes phase χ dependent and is given by

|G S,deg | 2 (χ) = 1+4 r 2 [ r -δ r cos(2χ)] + (1 + 2 r -δ 2 r ) sin(2χ) (1 -2 r + δ 2 r ) 2
.

(6) Equations (4-6) are valid only below the onset of parametric oscillations, that is of pump-induced autooscillations at zero signal B S for r > 1+δ 2 r . For sufficiently large pumping strength r > 0.42 the power gain |G S | 2 is larger than 2 at small ∆ r and δ r , and a gain bandwidth ∆ω = 2π∆f at -3 dB can be defined. For the optimal pumping frequency δ r = 0 we find

∆ω 2Γ a = (1 + 2 r ) 2 r 4 r -6 2 r -1 -1 , (7)
which yields a gain bandwidth product

|G S | ∆ω/2Γ a 1 (8)
that is constant within 10% above 7 dB gain. Then, saturation can be evaluated approximately in an easy way by noticing that as the internal amplitudes of oscillation A S and A I increase with the pumping strength and gain, they tend to the same value when

|G S | 2 |G I | 2
1 (see Eq. ( 5)). Consequently, the terms in α r in the first two equations of system (3) also converge to close values and play the very same role as the pump to resonator detuning δ r , which is itself responsible for a gain drop given by Eq. ( 6). Equating α r A 2 S + 2A 2 I at δ r = 0 to the value δ r,sat = 0.35 1 -2 r that produces a -1 dB drop of |G S | 2 -1, leads to the following equivalent values for A S and B S (so called 1 dB compression point):

   A 2 S,sat 0.35 |α R | √ 1-2 r 1+2 2 r B 2 S,sat Γa 0.17 |α R | (1-2 r ) 5/2 1+2 2 r (9)
In addition, saturation at large gain |G S | has to occur when the peak current i in the junctions is still well below their critical current i c . Since at δ r ∼ 0,

i/i c = (1 + r ) |A S | p √ Z r /N , keeping i/i c < 0.5 yields the design rule pQ > 21/ |G S,max | + 1, (10)
which imposes a minimum p for low Q and wide bandwidth JPAs.

In this work, we choose to implement a tunable amplifier in the 5-6 GHz range with a quality factor Q of order 100, which should have a product gain×bandwidth of ∼ 50 MHz according to Eq. ( 8). To reduce the maximum microwave pumping power corresponding to r = Qa P /2 = 1, i.e. to a modulation a P ∼ 1% of the total inductance, a high participation ratio p 0.5 is chosen. On the other hand, in order to keep the nonlinearity α r weak and to increase the 1 dB compression point B S,sat , the tunable inductance is implemented with N = 8 SQUIDs. In this case, the saturation power is increased by N 2 or 18 dB, compared to the case of a single SQUID with the same total Josephson inductance. Finally, as Z r plays only a minor role in the nonlinearityα r (in comparison with p 3 and N -2 ), its value will be simply chosen at the best convenience for implementing the lumped element resonator. The device was fabricated on a thermally oxidized Si chip by sputtering 170 nm of niobium and patterning the whole structure (except the SQUID array) by optical lithography and CF 4 -Ar reactive ion etching. The SQUID array was then fabricated by e-beam lithogra-phy and double-angle evaporation of aluminum with oxidation of the first Al layer. Each SQUID has a loop area of 8 µm × 15 µm and two junctions with nominal area 2.2 µm × 0.7 µm and tunnel resistance 141 Ω, yielding β 0.1. The active antenna wires of the pump line are positioned 16 µm away from the SQUID centers.

The measurement set-up is schematized in Fig. 1b. A small superconducting coil is used to compensate the global DC flux offset. DC flux biasing and AC pumping of the SQUIDs are obtained by two attenuated and filtered lines combined with a bias-tee. The input line includes attenuators at various temperatures and a 4 -8 GHz bandpass filter. The -71.5 dB and -51 dB transmissions of the input and pumping lines are calibrated with a ±1 dB uncertainty. The reflected and amplified signal is routed to the output line by a cryogenic circulator with -18 dB isolation. This output line includes an isolator for protecting the sample from higher temperature noise, a 4 -8 GHz filter, a cryogenic high electron mobility transistor (HEMT) amplifier at 4K with 38 dB gain and a calibrated noise temperature of 3.8 K, as well as additional room temperature amplifiers. The output signal is finally analyzed using a spectrum analyzer or a homodyne demodulator followed by a digitizer. Microwave generators for the input signal and pump are precisely phase locked.

III. EXPERIMENTAL RESULTS

Measurements were performed in a dilution refrigerator at the temperature of 30 -40 mK. As a preliminary characterization, the resonance was measured with a vectorial network analyzer by recording the phase of the reflected signal at zero pumping and at nominal input power P S,n = -126 dBm small enough to avoid any nonlinear effects (all nominal powers mentioned here and below refer to powers at the sample ports given the calibration of the lines). Inset of Fig. 2 shows this resonance at zero flux with a fit of the curve yielding a maximum frequency of f R0 = 5.97 GHz. The comparison with the f R1 = 8.06 GHz resonance frequency of a similar resonator with shorted junctions yields p = 1-(f R0 /f R1 ) 2 = 0.45, close to the 0.42 design value. Fitting the expression -2 arctan[2Q(f S /f R -1)] to the measured resonance curve also gives the quality factor Q 0 70, with however limited accuracy due to a setup imperfection yielding spurious multiple wave interferences (see shoulders in inset of Fig. 2). The main graph of Fig. 2 shows the variation of f R as a function of the applied flux Φ and its comparison with the theoretical prediction from section I. The agreement is only qualitative especially above 0.35 Φ 0 where f R decreases faster than predicted by our simple model that does not include either the flux inhomogeneity in the different SQUIDs or the possible penetration of the flux through the junctions.

For characterizing amplification, the input and output lines are then connected as shown in Fig. 1. The signal and idler gains are measured with the spectrum analyzer by comparing the output powers of the signal and idler without and with parametric pumping at f P = 2f R . The gains increase with |Φ DC | and the slope of the modulation curve at fixed absolute pumping power.

In the following measurement, signal amplification is fully characterized at the working point (Φ 1 /Φ 0 = -0.32, f R1 = 5.17 GHz), i.e. point B on Fig. 2, where the slope κ 1 = 1.62 is at the same time large and in agreement with the predicted value. At this point the SQUID array model predicts a participation ratio p 1 = 0.59 and a quality factor Q 1 = 81.

The non-degenerate (∆ = 0) signal power gain |G S | 2 is measured with the pump frequency f P = 2f R (δ = 0) as a function of the signal frequency f S for increasing nominal pump power P P,n , at an input power P S,n = -142.5 dBm sufficiently low to avoid the saturation at the highest gain. Close to the resonance, a minimum detuning δ/2π = 5 kHz is used to avoid operation in the degenerate mode. Figure 3a shows the gain increase up to 31.8 dB (dashed top curve, for which parametric oscillations are about to start) and the corresponding bandwidth decrease. The maximum power gain |G S | 2 and the corresponding -3 dB bandwidth ∆f deduced from Fig. 3 are plotted on Fig. 4b together with the amplitude-gain × bandwidth product |G S | ∆f . This product happens to be almost constant around 61 MHz over the whole 7 dB -30 dB gain range. Besides, the idler gain (data not shown) approaches the signal gain at large values.

In order to check that the amplifier operates close to the quantum limit, i.e. with a noise temperature of order T N = hf R /2k B 125 mK [22], the variation of the signal and noise powers are compared when switching on and off the pump: from the 2.9 dB increase of the noise when switching on a 18.4 dB gain, from the calibrated 3.8 ± 0.3 K noise temperature of the HEMT amplifier alone in a separate run, and from the 1.7 ± 0.2 dB attenuation of elements placed below 250 mK between the sample and the HEMT amplifier, we deduce an apparent noise temperature of only 80 ± 10 mK. This value is smaller than the expected quantum limit of 125 mK, a discrepancy that shows that modeling the line by a simple attenuator is not sufficient, as supported by our observation of the setup imperfection already mentioned. This result nevertheless indicates that our JPA is not far from the quantum limit. A more precise determination of T N would require a much more precise control and calibration of the low temperature part of the measurement line, as well as a switch to connect the detection chain either to the JPA or to a low temperature reference noise source [10].

The phase dependent gain in the degenerate case (δ = 0 ) was then measured with f S = 5.17 GHz by varying the phase χ of the signal with respect to the pump; it is shown on Fig. 3b for the two values of the pump power that correspond to a 20 dB gain and to the maximum gain in the non-degenerate mode. As expected the maximum degenerate gain is 6 dB larger than the nondegenerate gain at almost the same frequency. At the highest degenerate gain of 36.9 dB, it was checked using the IQ demodulator (see inset of Fig. 3b) that the phase of the amplified signal is stable over minutes and that the output signal drops down to zero (no parametric oscillation) when the input signal is switched off. As the phase χ is varied, the measured degenerate gain varies as expected, the lowest value of -25 dB resulting from the uncontrolled interference between the deamplified signal and the ∼ -18 dB leak of input signal through the circulator (see Fig. 1). This strong deamplification and the low noise temperature indicate that our JPA could also be used as a vacuum squeezer. In the inset of Fig. 3 the elongation of the Gaussian spot along the amplified quadrature shows that after parametric amplification the noise coming from the sample at f S = 5.17 GHz overcomes the noise of the cryogenic amplifier placed at 4 K, the size of which is given by the spot size in the perpendicular direction. In this latter direction we observe that the spot size is reduced by 1.1% when switching on the parametric pumping. This reduction is twice as small as the 2.2% expected from deamplification of vacuum noise, which is again related to the difficulty to determine the noise temperature of the whole setup.

Finally, the saturation of the JPA is measured by recording the non-degenerate signal power gain as a function of the signal input power P S,n for the same series of pump powers P P,n as before (see Fig. 4a). The signal gain is almost constant at low input power and then decreases above a P P,n dependent threshold in P S,n (however with a small bump of up to 1 dB just before saturation possibly due to higher orders in non-linearity). In practice, the 1 dB compression point is defined as the input power P S,sat at which the gain is 1 dB lower than at P S,n = -150 dBm; it is plotted on Fig. 4c. The set of measurements of Figs. 4b-c is then compared to the linear model of section I: the power gain, bandwidth, product amplitude-gain × bandwidth, as well as the 1 dB compression point of Eqs. (6)(7)(8)(9) are calculated by using the values of f R1 , p 1 , κ 1 and Q 1 indicated above and are plotted on Fig. 4. Given the ±1 dB uncertainty on the calibration of the signal and pump lines, the nominal input and pump powers were shifted by +0.8 dB and -1.0 dB to match the theory at the lowest pumping power. The agreement between the overall measured data and the model is surprisingly good given the crudeness of the linear model. This fair agreement validates the idea of increasing the number N of SQUIDs to increase the saturation power that scales with N 2 . With a single SQUID having the same total inductance as the array implemented here (about 1.7 nH), the saturation would have been N 2 = 18 dB lower, as indicated by the dashed line of Fig. 4c. The discrepancy between experimental data and the model increases with P p as the nonlinearity plays a more important role, and the actual parametric amplification region extends a bit over the theoretical parametric oscillation region of the linear model (dotted line of Fig. 4). The performances of the present device are comparable to those of other truly parametric amplifiers recently made. Due to our choice of a rather large Q ∼ 70, the gain bandwidth product is smaller than what was obtained for instance in [23] with Q ∼ 10. In [15], the direct coupling of the resonator to a cleverly engineered frequency dependent external impedance yielded an even lower Q and a bandwidth above 500 MHz. Despite the use of N = 8 SQUIDs, the 1dB compression point obtained here is not very high due to its scaling as N 2 Q -2 p -3 and to the large participation ratio and quality factors chosen to minimize the pump power: it is however about 12 dB above a similar amplifier made of a single SQUID with about the same critical current [18] and only a few dB below another one [23] with smaller participation ratio p (three times larger critical current Ic) and Q.

In terms of perspectives, equations (8-10) predict that with a similar geometry N∼ 10, a smaller Q ∼ 10, and higher critical currents yielding p ∼ 0.25, a bandwidth of ∼ 50 MHz, and a compression point ∼ -100 dBm should be obtained at 20 dB gain. This would require a larger pump power, i.e. a larger flux modulation Φ AC ∝ 1/pQ at constant gain, which would reach 0.1Φ 0 . Such a large modulation could be technically difficult to achieve. Increasing the number of SQUIDs is also an obvious optimization axis: If theoretically the array length has just to be kept much smaller than the pump wavelength so that all SQUIDS are pumped in phase, the practical difficulty is to DC flux bias and homogeneously modulate all the SQUIDs.

In summary, a lumped element truly parametric Josephson amplifier has been designed and characterized. Its inductance is implemented by a SQUID array to limit its nonlinearity and increase the maximum allowed input power. With a quality factor of 70 -80, this simple device provides a gain of up to 30 dB, a product amplitude-gain×bandwidth of 61 MHz, and a 1 dB compression point of -123 dBm at 20 dB gain. Although its behavior is in agreement with theory and demonstrates the advantage of using a SQUID array, it can still be optimized by reducing both its quality factor and its Josephson participation ratio to the inductance and/or by increasing the number of SQUIDs in the array. Operated close to the quantum limit, this truly parametric amplifier could also be used as a quiet and strong squeezer in degenerate mode or as the first stage of amplification in a superconducting quantum bit readout.

  II. SAMPLE AND MEASUREMENT SETUPAn optical micrograph of the JPA and its equivalent circuit are shown in Fig.1a. This JPA is made of an interdigitated coplanar capacitor to ground (split in two parts) with capacitance C R = 0.40 pF , in parallel with an inductance to ground L combining in series a meander of inductance L g = 0.80 nH with an array of 8 SQUIDs of total Josephson inductanceL J = 0.88 nH at zero magnetic flux Φ. Being designed to be operated in reflection, this LC circuit is coupled to a single inputoutput signal line (50 Ω coplanar waveguide -CPW ) through a C c = 55 fF capacitance yielding a characteristic impedance Z r 65 Ω and a quality factor Q 70 at Φ = 0. On the other side of the device a 50 Ω CPW line shorted to ground by two loops coupled inductively to 4 SQUIDS each, serves both for their DC flux biasing and for parametric pumping. Note that after compensation of any global DC flux offset, the magnetic fluxes Φ are exactly opposite in the left and right 4 SQUID sub-arrays, which yields the same inductance modulation.

Figure 1 :

 1 Figure 1: Experimental setup. (a) Optical micrograph of the tested parametric amplifier showing its 50 Ω coplanar waveguide (CPW) signal input port (top), its coupling capacitor Cc, its capacitor CR (left and right), its inductor Lg (middle) terminated by an 8-SQUID arrays with total inductance LA (Φ), and its split magnetic flux line coupled to a 50 Ω CPW (bottom). The DC current in the flux line sets the dc flux ΦDC and the resonance frequency fR of the resonator, whereas the ac current parametrically pumps the resonator at ωP 2ωR. The black arrow points to the equivalent circuit. (b) Electrical circuit diagram showing (from left to right) the pump line, the dc flux line added to the pump with a bias tee, an additional flux line feeding a coil for compensating any flux offsets, the signal line, a circulator routing the reflected and amplified signal to a measurement line through an isolator protecting the sample from the noise of the first amplifier placed at 4 K. Feeding lines are attenuated and filtered. The output signal is split after amplification and analyzed both with a spectrum analyzer and by homodyne demodulation.

Figure 2 :

 2 Figure 2: DC flux modulation. Experimental (dots) and calculated (line) resonator frequency fR as a function of the DC current in the on-chip flux line (see Fig.1a) measured by fitting the phase of a weak signal reflected on the resonator in absence of parametric pumping, as shown in the inset for zero flux bias (point A). We attribute the two shoulders on the sides of the measured resonance to multiple wave interferences due to an imperfect impedance matching somewhere in the setup. Parameters used for calculation of the modulation curve are fR0, p = 0.45, and β = 0.1 . Full characterization in the next figures are done at working point B.

Figure 3 :

 3 Figure 3: Signal power gain |GS| 2 at working point B of Fig. 2 for a nominal input power PS,n = -142.5 dBm. (left) Nondegenerate gain as a function of the signal frequency fS at different nominal pumping powers PP,n between -58.7 dBm and -54.4 dBm (top dashed curve just before the onset of parametric oscillation in the absence of incident signal). (b) Phase sensitive degenerate gain for PP,n = -54.4 dBm and -55.4 dBm. (inset) Demodulated signal in the IQ plane at maximum gain (PP,n = -54.4 dBm) filtered at 1 MHz. I and Q voltages are digitized at 1 MSample/s during 2 s, and the color encodes the density of samples from 0 (dark blue) to maximum (red).

Figure 4 :

 4 Figure 4: Amplifier characterization at the working point B of Fig. 2 for non-degenerate pumping. (a) Signal power gain as a function of the nominal input power PS,n showing the saturation at the same nominal pumping powers PP,n as in Fig. 3a (top dashed line corresponds again to the onset of parametric oscillation). Vertical dashed and dotted line correspond to the input powers where gain was measured in Fig. 3 and where reference gain for saturation was defined, respectively. (b) Power gain |GS| 2 , bandwidth BW , and product |GS| × BW deduced from measurements (dots) of Fig. 3a at PS,n = -142.0 dBm, and calculated (solid lines) from the model with the parameters indicated in the text. (c) 1 dB compression point deduced from (a) (dots), calculated from the model (solid line), and calculated with the sames parameters but only one SQUID (dashed line). Note that given the ±1 dB precision on the calibration of the pumping and signal lines, +0.8 dB and -1.0 dB were added to the nominal PP,n and PS,n values to match the data to the theoretical curves at low pumping strength. The vertical dotted line indicates the frontier between parametric amplification and parametric oscillation (infinite gain) for the linear model.

Acknowledgment

We gratefully acknowledge discussions within the Quantronics group, technical support from P. Orfila, P. Senat, J.C. Tack and Dominique Duet, as well as financial support from the European research contract SCALE-QIT.


Carles Altimiras

Olivier Parlavecchio

Philippe Joyez

Denis Vion

Patrice Roche

Daniel Esteve

Fabien Portier

Dynamical Coulomb Blockade of Shot Noise

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

Dynamical Coulomb Blockade of Shot Noise

Carles Altimiras, * Olivier Parlavecchio, Philippe Joyez, Denis Vion, Patrice Roche, Daniel Esteve, and Fabien Portier † Service de Physique de l'Etat Condensé (CNRS URA 2464), IRAMIS, CEA-Saclay, 91191 Gif-sur-Yvette, France (Dated: September 5, 2013)

We observe the suppression of the finite frequency shot-noise produced by a voltage biased tunnel junction due to its interaction with a single electromagnetic mode of high impedance. The tunnel junction is embedded in a λ/4 resonator containing a dense SQUID array providing it with a characteristic impedance in the kΩ range and a resonant frequency tunable in the 4-6 GHz range. Such high impedance gives rise to a sizeable Coulomb blockade on the tunnel junction (∼30% reduction in the differential conductance) and allows an efficient measurement of the spectral density of the current fluctuations at the resonator frequency. The observed blockade of shot-noise is found in agreement with an extension of the dynamical Coulomb blockade theory.

Contrarily to usual electronic components for which one can define an intrinsic behavior (e.g. the I-V characteristic), the transport properties of a coherent quantum conductor depend on its biasing circuit. This is true even when the size of the circuit exceeds the electron coherence length, suppressing electronic interference effects. This non-intrinsic behavior can be traced to the quantumprobabilistic character of the transmission of electrons through the conductor, resulting in broad-band fluctuations of the current called shot noise [? ]. This current noise can create collective excitations (hereafter called "photons") in the electromagnetic environment seen by the conductor. This yields a back action on the transport properties of the conductor itself [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF]. This physics bears similarities with the spontaneous emission of photons by an excited atom, albeit with important differences: first, dc biased quantum conductors are out of equilibrium open systems and cannot be described as a set of discrete levels. Second, the dimensionless parameter characterizing the electron-photon coupling is given by the ratio of the environment's impedance to the resistance quantum R K = h/e 2 25.8kΩ; Hence, by increasing the impedance of the circuit connected to the quantum conductor, one can significantly increase the effective coupling constant. This results in a rich physics, already partially understood: noticeably, the dynamical Coulomb Blockade (DCB) theory [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF] accounts for the observed suppression [3][4][5] of the low voltage conductance of a tunnel element as a result of its coupling to a dissipative electromagnetic environment. A natural step is then to understand how the coupling to the environment modifies the current fluctuations themselves: is there a Coulomb Blockade of shot noise? This question of current fluctuations in the presence of DCB was addressed theoretically for the low frequency/long time limit, where the corrections to the noise power and to the full counting statistics were predicted [START_REF] Galaktionov | Current fluctuations and electron-electron interactions in coherent conductors[END_REF][7][8][9]. Instead, we consider here the environment feedback on the frequency dependence of the shot noise of a simple quantum conductor, a tunnel junction. Following Ref. [10], we extend the DCB theory to predict the finite frequency emission noise spec-trum of a voltage biased tunnel junction in the presence of an arbitrary linear environment. Probing this prediction requires to achieve strong coupling of the junction to its environment and to measure its high frequency shot noise. To do so, we fabricate a tunnel junction embedded in the simplest environment, a harmonic oscillator, and measure the effect of Coulomb blockade on the shot noise power at the frequency of the oscillator. The oscillator is realized with a microwave resonator based on a Josephson transmission line allowing both a tenfold increase of the coupling constant between the junction and the resonator, and to tune the resonant frequency. The data are found in quantitative agreement with the theory.

In order to evaluate the current and its fluctuations through a tunnel element in the presence of DCB, we consider a circuit consisting of a tunnel junction of conductance G T in series with an impedance Z(ν) described as the series combination of harmonic modes (see upper panel of Fig. 1) at temperature T and biased at voltage V . We then compute (see the Supplemental Material for more details) the current I and the quantum spectral density S I (ν) of current noise, i.e. the Fourier transform of the non-symmetrized current-current correlator

S I (ν) = 2 ∞ -∞ I(t)I(0) e -i2πνt dt. (1)
In this convention positive (resp. negative) frequencies correspond to energy being emitted (resp. absorbed) by the quasiparticles to (resp. from) the electromagnetic modes. Taking separate thermal equilibrium averages over the unperturbed quasiparticle and environmental degrees of freedom yields

I(V ) = G T e [γ * P (eV ) -γ * P (-eV )], S I (ν, V ) = 2G T [γ * P (eV -hν) + γ * P (-hν -eV )],
where γ * P (E) = dε γ(ε )P (E -) with P (ε) the probability density for a tunneling electron to emit the energy ε in form of photons into the impedance [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF], with γ(

) = dε f (ε )[1 -f (ε + ε)] = ε/(1 -e -ε/k B T
), and with f the Fermi function. Eq. 2 is the standard DCB expression for the tunneling current [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF], whereas Eq. 3 is our prediction for the Coulomb Blockade of shot noise, which we probe in the experiment described below. For a positive bias voltage and low temperature (k B T eV, hν), Eqs. 2-3 take the simpler form

I(V ) = G T e eV 0 (eV -ε)P (ε)dε, (2)
S I (ν, V ) = 2G T Θ(eV -hν) eV -hν 0 (eV -hν -ε)P (ε)dε, ( 3
)
with Θ(ε) the Heavyside function. Eqs. 4-5 are easily interpreted: the total current is proportional to the average energy available for quasiparticles upon the transfer of an electron through the circuit, and so is the noise power at frequency ν, albeit imposing the emission of a photon of energy hν into the environment. In the case of a vanishing impedance Z(ν) R K , P (ε) = δ(ε) and one recovers the standard, non interacting, finite frequency shot noise result [START_REF] Lesovik | Here the distinction is made between emission noise SI(ν) = 2 ∞ -∞ I(t)I(0) e -i2πνt dt and absorbtion noise SI(-ν). While observation of the later requires excitation of the sample by external sources, for a zero temperature external circuit, only SI(ν) is observed. For an earlier high-frequency shot noise derivation not making the distinction between SI(ν) and SI(-ν)[END_REF]. In the case of a discrete harmonic oscillator of frequency ν 0 = 1/ 2π √ LC and impedance Z C = L/C, thermalized at a temperature T hν 0 /k B :

P (E) = ∞ k=0 p k δ(E -khν 0 )
, with p k = e -α α k /k! the probability for the oscillator to absorb k photons [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF], and α = πZ C /R K the coupling strength between the tunnel junction and the oscillator. Our experiment achieves an unprecedented electron-single mode coupling α ∼ 0.3, which allows observing multi-photon processes both in the average current and in the emission noise. Note that, despite similar denominations, the effect we consider here differs from static Coulomb Blockade, which results from the charging energy of a small island connected to reservoirs by tunnel barrier. Static Coulomb Blockade is a quasiclassical effect which can be described by master rate equations, at the level of the current noise [START_REF] Kafanov | Measurement of the shot noise in a single-electron transistor[END_REF], and even the full counting statitsics [START_REF] Gustavsson | Counting statistics and super-Poissonian noise in a quantum dot: Time-resolved measurements of electron transport[END_REF].

Our experimental set-up is schematized in the lower panel of Fig. 1: a 100 × 100 nm 2 tunnel junction with tunnel resistance G -1

T in the 100 kΩ range is embedded in an on-chip λ/4 coplanar resonator of resonant frequencyν 0 , whose inner conductor is made of an array of identical and equally spaced Al/AlOx/Al SQUIDs. To a very good approximation, its lineic inductance is dominated by the Josephson inductance L J = (2eI 0 | cos(eφ/ )|a) -1 , where I 0 is the maximum critical current of one SQUID, φ the flux applied to each SQUID, and a the distance between adjacent SQUIDs. This increases the resonator impedance Z C above 1 kΩ, and allows to decrease ν 0 while increasing Z C by applying a flux through the SQUIDs. Two samples were fabricated and measured. In both cases the 6 GHz maximum frequency of the resonator ensures k B T hν 0 at the 15 mK temperature of the experiment, so that thermal fluctuations do not blur Coulomb Blockade effects. The minimum zero flux lineic inductances of the first and second resonators were designed at 8.80 10 -5 Hm -1 and 3.95 10 -4 Hm -1 , respectively. Note that Josephson transmission lines have been used to create non-linear resonators used as parametric amplifiers [START_REF] Castellanos-Beltran | Amplification and squeezing of quantum noise with a tunable Josephson metamaterial[END_REF], or to probe how quantum phase slips drive them into an insulating state at Z C R K [START_REF] Chow | Length-Scale Dependence of the Superconductor-to-Insulator Quantum Phase Transition in One Dimension[END_REF]. We avoid this regime by keeping Z C in the kΩ range, which stills allows us to obtain sizeable DCB corrections. Keeping the current going through the resonator much smaller than I 0 ensures that the Josephson junctions can be considered as linear inductances. The resistance G -1

T being much higher than Z(ν), the impedance seen by one conduction channel of the junction is not shunted by the parallel conductance of the other channels [START_REF] Joyez | How Is the Coulomb Blockade Suppressed in High-Conductance Tunnel Junctions?[END_REF]. The SQUIDs and the tunnel junction were fabricated on a Si/Si0 2 substrate using standard nanofabrication techniques [START_REF]For more information[END_REF]. In addition, a 30 × 50 × 0.3 µm 3 gold patch is inserted between the tunnel junction and the SQUID array in order to evacuate the Joule power dissipated at the tunnel junction via electron-phonon coupling. As an example, assuming a typical 2nW µm -3 K -5 electron-phonon coupling constant [START_REF] Huard | Electron heating in metallic resistors at sub-Kelvin temperature[END_REF], a 100µV (resp. 1mV) bias on a 200kΩ tunnel resistance increases the electron temperature from 15mK to 20mK (resp 50mK), keeping heating effects negligible. Note that the thermalization pad adds an additional 12fF to ground, which is taken into account to evaluate the total impedance seen by the tunnel junction. The chip is connected to the biasing and measurement circuits through a commercial 50 Ω matched bias tee. The inductive (low frequency) path is used both to bias the sample through a cold 13 MΩ resistor, and to measure the dc voltage across the tunnel junction and its conductance G(V ). The capacitive (RF) path guides the radiation S I (ν, V ) emitted by the sample to a cryogenic isolator anchored at 15 mK, to a cryogenic amplifier with a ∼2.5 K noise temperature in the 4-8 GHz bandwidth, to room temperature band-pass filters, and finally to a power "square law" detector, the output voltage of which is proportional to its input microwave power. The isolator diverts the current noise of the amplifier to a 50Ω matched resistor that re-emits to the sample a blackbody radiation only at the coldest temperature, ensuring a negligible photon occupation of the resonator at GHz frequencies. Finally, the signal S I (ν, V ) is extracted from the large noise floor of the cryogenic amplifier by a lockin detection involving a 1 µV sinusoidal modulation at 17 Hz on top of the dc voltage V .

We first characterized the on-chip microwave resonator by measuring the power emitted by the electronic shot noise of the junction S I ∼ 2eI at high bias voltage V ∼ 1 mV [START_REF] Hofheinz | The Bright Side of Coulomb Blockade[END_REF], where DCB effects are negligible. Under these conditions, the spectral density of the emitted power is

2eV Re[Z(ν)]G T /|1 + G T Z(ν)| 2 2eV G T Re[Z(ν)] since the tunnel resistance (G -1
T =230 kΩ /450 kΩ for sample 1/2) is much larger than the maximum detection impedance Z(ν). This spectral density is obtained using a heterodyne measurement implementing a 10 MHz-wide band pass filter at tunable frequency. As shown in Fig. 2, the extracted Re[Z(ν)] is in satisfactory agreement with predictions. In particular, Z(ν) shows the expected resonance, with a resonant frequency ν 0 decreasing with φ, associated to an increasing impedance and quality factor, which is limited by radiative losses. The maximum disagreement between the measured maximum for Re[Z(ν)] and the calculated one is about 15%, which we attribute to an uncertainty in the calibration of the gain of the amplifying chain [START_REF]For more information[END_REF]. We attribute the additional structure around 5.7 GHz to a parasitic resonance in the detection chain. Once our microwave environment calibrated, we measure both the differential conductance G(V ) of the tunnel junction and the voltage derivative ∂S I (ν 0 , V )/∂V of the noise emitted in a 180 MHz bandwidth centered around the resonator frequency ν 0 , as a function of the dc bias voltage applied to the junction. The conductance, shown in the upper panel of Fig. 3 for both samples, is non-linear, showing a stair-case behavior characteristic of DCB corrections due to a single mode, rounded by the finite temperature [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF]. The high characteristic impedance of our resonators yields DCB corrections to the conductance 10 times higher than with standard microwave resonators [START_REF] Hofheinz | The Bright Side of Coulomb Blockade[END_REF][START_REF] Holst | Effect of a transmission line resonator on a small capacitance tunnel junction[END_REF], and shows not only the single photon emission onset at bias voltage V 0 = hν 0 /e but also the two photon onset at 2V 0 . As shown in the lower panels of Fig. 3, the voltage derivative of emission noise power also displays a non-linear staircase shape, with a first singularity at V 0 , followed by a smaller step at 2V 0 . The first step at V 0 is predicted in the standard -non including DCB effects-finite frequency shot noise theory, represented by the dotted black curve in the lower panel of Fig. 3. It has been observed in several experiments [START_REF] Schoelkopf | Frequency dependence of shot noise in a diffusive mesoscopic conductor[END_REF][START_REF] Onac | Using a quantum dot as a high-frequency shot noise detector[END_REF][START_REF] Zakka-Bajjani | Experimental test of the high-frequency quantum shot noise theory in a quantum point contact[END_REF][START_REF] Gustavsson | Frequency-selective single-photon detection using a double quantum dot[END_REF] and can be dually understood either in

1 . 0 -2 -1 0 1 2 0 6 1 2 1 8 0 1 2 0 6 1 2 1 8 0 1 2 -3 -2 -1 0 1 2 3 -5 0 -3 3 -1 7 0 1 7 3 3 5 0 0 6 1 2 1 8 0 1 2 Sample 1: G T =3.72 µS G / G T V (µV) ( 2 e G T ) -1
. terms of the finite time coherence of a DC biased quantum conductor, or in terms of the energy cost of creating excitations at frequency ν 0 in the measuring apparatus [START_REF] Lesovik | Here the distinction is made between emission noise SI(ν) = 2 ∞ -∞ I(t)I(0) e -i2πνt dt and absorbtion noise SI(-ν). While observation of the later requires excitation of the sample by external sources, for a zero temperature external circuit, only SI(ν) is observed. For an earlier high-frequency shot noise derivation not making the distinction between SI(ν) and SI(-ν)[END_REF][START_REF] Gavish | Detection of quantum noise[END_REF][START_REF] Aguado | Double quan-tum dots as detectors of high-frequency quantum noise in mesoscopic conductors[END_REF]. The second step occurs at the onset voltage for the emission of two photons in the resonator by a tunneling electron. The significant difference between the experimental points and the non interacting prediction demonstrates the Coulomb Blockade of shot noise.

ٖS I / ٖV eV/hν 0 ν 0 = 6GHz ν 0 = 6GHz ν (GHz) Sample 2: G T = 2.22 µS Z C ( k Ω) ν (GHz) Z C ( k Ω) ν 0 = 4GHz Z C ( k Ω) ν (GHz)
We now probe how the data shown in Fig. 3 can be quantitatively accounted for by Eqs. 2 and 3, using our well controlled environment as an input to evaluate P (E). We model this environment as a series combination of three discrete harmonic modes. The two higher frequency ones correspond to the fundamental and first harmonic modes of the resonator. These two modes account with no adjustable parameters for the observed variations above V 0 = hν 0 /e. Their characteristic impedance can be evaluated through the standard formula Z C = 2 ν0ImY (ν0) , where Y (ν) is the environment's admittance, evaluated from our modeling of the Josephson transmission line. We introduce an additional lower frequency mode to account for the unexpected 3% dip in the differential conductance that we observe at low bias voltage |V | 5µV . We attribute this low-frequency parasitic resonance, which only slightly affects the data, to the bias-Tee. The corresponding predictions, assuming an electron temperature T e = 16 mK corresponding to the temperature of the fridge's mixing chamber, are represented by the solid red curve in the top graphs of Fig. 3. Note that at this temperature, the ∼ 3.5k B T /h ∼ 1 GHz smearing expected from the Fermi distribution is broader than the linewidth of the modes of our resonator. This is why the discrete modes model, which yields analytical expression for the P (E) function [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF], is able to reproduce the data. The emission noise data can be reproduced by Eq. 2 with excellent accuracy, whereas the expression corresponding to the current noise spectral density symetrized with respect to frequency [10], S sym I (ν, V ) = [S I (-ν, V ) + S I (ν, V )]/2, represented by the green dashed-dotted line in Fig. 3, is not compatible with our data [START_REF]For more information[END_REF]. Note that at low tem-perature k B T hν 0 the relative size of the two-photon step is α/2, which explains why noise blockade is not seen with usual environment impedance yielding values of α ∼ 10 -3 . However, when considering the proposed primary shot noise thermometry [START_REF] Spietz | Primary Electronic Thermometry Using the Shot Noise of a Tunnel Junction[END_REF], even such low values of α cause a systematic correction that should be considered to reach metrological accuracy. Last, the data shown in the to left panels of Fig. 3 were taken for the maximum value of ν 0 = 6 GHz. Applying a flux through the SQUIDs induces a stronger blockade due to the increased Josephson inductance. As shown on the right panel of Fig. 3, the higher impedance of sample 2 allows observing even three-photon processes when the resonant frequency is set at 4 GHz (the lower end of our detection bandwidth), yielding Z C 2.25 kΩ and a 30% reduction of the zero-bias conductance.

In conclusion, we have developed an original electromagnetic environment, allowing to reach an unprecedented coupling between a quantum conductor and a single mode environment. We took advantage of this to demonstrate the Coulomb Blockade of the finite frequency noise of a tunnel junction. Two-and threephoton processes are identified, in agreement with an extension to the theory of Dynamical Coulomb Blockade. The experimental methods developed here can be readily applied to quantum conductors of arbitrary transmissions, for which a complete description of quantum transport in the presence of an electromagnetic environment is still missing. Noticeably, they allow to probe the Coulomb Blockade of shot noise in quantum point contacts [START_REF] Golubev | Coulomb interaction and quantum transport through a coherent scatterer[END_REF][START_REF] Levy Yeyati | Direct link between coulomb blockade and shot noise in a quantum-coherent structure[END_REF][START_REF] Cron | Electronic Correlations: From Meso to Nano-Physics[END_REF][START_REF] Altimiras | Experimental test of the dynamical coulomb blockade theory for short coherent conductors[END_REF][START_REF] Souquet | Dynamical Coulomb blockade in an interacting one-dimensional system coupled to an arbitrary environment[END_REF][START_REF] Parmentier | Strong back-action of a linear circuit on a single electronic quantum channel[END_REF], where DCB was recently demonstrated to bear a deep connection to the physics of impurities in Luttinger liquids [START_REF] Jezouin | Tomonaga-Luttinger physics in electronic quantum circuits[END_REF], or quantum dots, where the interplay between resonant tunneling through the dot and the coupling to the environment was mapped to the physics of Majorana fermions [START_REF] Mebrahtu | Quantum phase transition in a resonant level coupled to interacting leads[END_REF]. This project was funded by the CNano-IDF Shot-E-Phot and Masquel, the Triangle de la Physique DyCoBloS and ANR AnPhoTeQ grants.

Technical assistance from Patrice Jacques, Pierre-François Orfila and Pascal Sénat, as well as discussions within the Quantronics group, with Inès Safi, Pascal Simon and Jean-René Souquet are gratefully acknowledged. * Presently at NEST, Istituto Nanoscienze-CNR and Scuola Normale Superiore, I-56127 Pisa, Italy † fabien.portier@cea.fr [START_REF] Ingold | Charge tunneling rates in ultrasmall junctions[END_REF]

CALCULATION OF THE CURRENT NOISE

The circuit we consider consists of a pure tunnel element connected in series with an arbitrary linear electromagnetic environment (See Fig. 1 of the main text). In this picture, the capacitance of the real tunnel junction has been incorporated into the impedance of the electromagnetic environment [START_REF] Ingold | Charge tunneling rates in ultrasmall junctions[END_REF]. The Hamiltonian of the circuit is:

H = H 0 + H T .
with

H 0 = ε c † c + r ε r c † r c r + H env .
Here, the indexes and r span all quasiparticle states in the left and right electrodes, the c † ,r (resp. c ,r ) denote the fermionic quasiparticle creation (resp. destruction) operators, H env is the Hamiltonian of the electromagnetic environment, and H T = T +T † is the tunneling Hamiltonian with T and T † implementing the transfer of one electron across the barrier from left to right and from right to left, respectiveley. This operator can be decomposed as T = e iϕ Θ where the e iϕ operator acts only on the electromagnetic environment by translating the charge transfered through the impedance by e, the charge of the tunneling electron, while Θ = ,r τ r c † c r acts only on the quasiparticles in the electrodes. In this writing τ r is the tunnel coupling amplitude of states and r, and ϕ denotes the phase difference operator across the tunnel element, acting on the electromagnetic environment, and related to the voltage drop V across the tunnel element by V = e ∂ϕ ∂t . The current operator through the tunnel element is given by I = -e ∂H ∂ϕ = -i e (T -T † ). The non-symmetrized noise current density S I (ν) is the Fourier transform of the current-current correlator:

S I (ν) = 2 ∞ -∞ I(t)I(0) e -i2πνt dt (1)
In this convention positive (resp. negative) frequencies correspond to energy being emitted (resp. absorbed) by the quasiparticles to (resp. from) the electromagnetic modes, and the current-current correlator function reads

I(t)I(0) = e 2 2 { T (t)T † (0) + T † (t)T (0) }.
We evaluate the averages by taking separate thermal equilibrium averages over the unperturbed quasiparticles and environmental degrees of freedom:

I(t)I(0) = e 2 2
e iϕ(t) e -iϕ(0) Θ(t)Θ † (0) + e -iϕ(t) e iϕ(0) Θ † (t)Θ(0) .

(

) 2
Unless at zero bias voltage where Eq. 2 is always exact, this amounts to assuming that there exist unspecified relaxation mechanisms fast enough to restore thermal equilibrium both in the electrodes and in the impedance between tunneling events. The validity of this latter assumption is highly dependent on the details of the system and should be checked on a case by case basis. Within this assumption, the occupation probability of the quasiparticle energy levels in the electrodes is given by the Fermi function f ( ) at the inverse temperature β. We further assume constant densities of states ρ l , ρ r in the electrodes and we replace the |τ r | 2 by their average value |τ | 2 over all l, r states. The average over quasiparticles degrees of freedom thus takes the form

Θ(t)Θ † (0) = ,r |τ r | 2 c † (t)c (0)c r (t)c † r (0) = |τ | 2 ρ l ρ r dεdε f (ε)(1 -f (ε + ε ))e -iε t/ = G T 2πe 2 dε γ(ε )e -iε t/ = G T 2πe 2 γ(t),
where G T = 4π 2 G K |τ | 2 ρ l ρ r denotes the tunnel conductance of the junction, which we assume much smaller than G K = e 2 /h, allowing to treat H T at the lowest order in perturbation theory; γ(t) denotes the inverse Fourier transform of

γ( ) = dε f (ε )(1 -f (ε + ε)) = 1-e -β = 2 1 + coth β 2 = n B (| |) + θ( ),
where n B is the Bose distribution function. The function γ( ) is proportional to the number of possible electronhole excitations in the Fermi seas of the electrodes with an energy difference of between the final and initial states. The same function is encountered in other contexts, such as in the noise of a resistor, for emission and arXiv:1403.5999v2 [cond-mat.mes-hall] 23 Oct 2014

absorption by bosonic degrees of freedom (see below). Even though γ( ) is not bounded, its inverse Fourier transform can nevertheless be expressed using distribution functions, meant to be integrated with proper functions [START_REF] Joyez | Single-electron tunneling at high temperature[END_REF]3]:

γ(t) = iπ 2 d dt δ(t) - π 2 β 2 sinh -2 πt β . (3)
The second quasiparticle term Θ † (t)Θ(0) can be checked to give the same result, i.e. Θ † (t)Θ(0) = Θ(t)Θ † (0) . To evaluate the average over the electromagnetic degrees of freedom, we decompose the phase difference across the tunnel element ϕ(t) = ϕ 0 (t) + φ(t) into the deterministic part ϕ 0 (t) = eV t/ corresponding to the dc voltage V across the junction, and a fluctuating random phase φ(t) caused by the fluctuations in the electromagnetic environment. For a linear electromagnetic environment at equilibrium, phase fluctuations are Gaussian. In this case, the averages we need to evaluate in Supplemental Material Eq. 2 can be expressed as exp ±i φ(t) exp ∓i φ(0) = exp J(t), where J(t) = ( φ(t)φ(0)) φ(0) is the phase-phase correlation function [START_REF] Ingold | Charge tunneling rates in ultrasmall junctions[END_REF]. The phase being proportional to the time derivative of the voltage, the phase-phase correlation function is related to the voltage noise on impedance seen by the junction, which can be obtained by the quantum fluctuation-dissipation theorem [START_REF] Ingold | Charge tunneling rates in ultrasmall junctions[END_REF]:

J(t) = hR -1 K +∞ -∞ S V (ν) (hν) 2 (e -i2πνt -1)dν, (4)
with

S V (ν) = 2ReZ(ν)γ(hν)
being the (equilibrium) non-symmetrized voltage noise across the environment impedance [4]. Collecting the above terms in Supplemental Material Eq. 2, we obtain the final expression for the current-current correlator:

I(t)I(0) = G T 2π 2 cos eV t e J(t) γ(t).
Inserting this result in the Fourier transform (1), the current noise can be expressed as

S I (ν) = 2G T [γ * P (-hν + eV ) + γ * P (-hν -eV )],
where * denotes convolution, and

P (E) = 1 2π ∞
-∞ e J(t)+iEt/ dt is the Fourier transform of exp[J(t)]. This function P (E) is interpreted as the probability for the environment to absorb the algebraical energy E during a tunnel event, and the convolution product:

γ * P (E) = dε γ(ε )P (E -)
simply accounts for all the possible ways to split an energy E between the quasiparticles degrees of freedom and the environment modes. For the case of an environment of vanishing impedance, P (E) → δ(E), and we recover a well known result

S I (ν) → S I (ν) = 2G T [γ(-hν + eV ) + γ(-hν -eV )].

SAMPLE FABRICATION

The 300 nm thick gold ground plane of the resonator and thermalization pad were obtained by optical lithography, followed by evaporation and lift-off. SQUIDs where fabricated following the process described in Ref. [5]: the SQUIDs (see the top inset) are obtained by double angle deposition of (20/40 nm) thin aluminum electrodes, with a 20 oxydation of the first electrode at 400 mBar of a (85% O 2 /15% Ar) mixture. Before the evaporation, the substrate was cleaned by rinsing in ethanol and Reactive Ion Etching in an oxygen plasma [START_REF] Galaktionov | Current fluctuations and electron-electron interactions in coherent conductors[END_REF]. The normal junction was obtained using the same technique, with 30/60 nm thick copper electrodes and an aluminum oxyde tunnel barrier (5 nm thick aluminum oxidized for 15 minutes at a 800 mBar (85%O 2 , 15%Ar) mixture).

DETAILS ON THE JOSEPHSON TRANSMISSION LINE

For sample 1 (resp. 2), our resonator consists in a 360 µm (resp. 91 µm) long Josephson metamaterial line containing 72 (resp. 38) lithographically identical and evenly spaced SQUIDs with a 5 µm (resp. 2.4 µm) period. The SQUIDs tunnel barriers of sample 1 (resp. sample 2) have an area of 0.5 µm 2 (resp. 0.5 µm 2 ) each resulting in a room temperature tunnel resistance R N = 720 Ω (resp. R N = 1880 Ω). To assess that the SQUIDs in the array are identical, we have performed reproducibility tests, yielding constant values of R N (within a few %) over millimetric distances. Assuming a superconducting gap ∆ = 180 µeV and a 17% increase of the tunnel resistance between room temperature and base temperature [7], one obtains a zero flux critical current for the SQUIDs I C = 671nA for sample 1 and I C =268 nA, corresponding to L J (φ = 0) = 0.49 nH for sample 1 and L J (φ = 0) = 1.25 nH for sample 2. This corresponds to an effective lineic inductance L 100 µH.m -1 at zero magnetic flux and frequency much lower than the Josephson plasma frequency of the junctions ν P [8]. Assuming a capacitance for the junctions of the order of 80 fF/µm 2 yields ν P 25 GHz. Note that our simple fabrication mask produces 10 times bigger Josephson junction in between adjacent SQUIDs, resulting in an additional ∼ L 10 µH.m -1 lineic inductance for sample 1 and ∼ L 25 µH.m -1 for sample 2. The ∼ L 1 µH.m -1 electromagnetic inductance associated to our geometry is negligible. With the designed lineic capacitance C = 75 pF.m -1 , the length of the resonator sets the first resonance at ν 0 8 GHz for sample 1 and ν 0

12 GHz for sample 2. The 12 fF shunting capacitance of the themalization pad reduces these frequencies to ν 0 6 GHz in both cases.

DETAILS ON SETUP AND CALIBRATION

We describe the calibration of the low frequency circuitry for voltage bias and current measurement as well the microwave components used to define the environment of the junction and to measure the emitted radiation.

Low frequency circuit

In addition to the components depicted in Fig. 1 of the main text, the low frequency circuit includes a copper powder filter anchored on the mixing chamber, as well as a distributed RC filter made with a resistive wire (50 cm of IsaOhm 304 Ω m -1 ) winded around a copper rod, and glued with silver epoxy on a copper plate in good thermal contact with the mixing chamber. Both are inserted between the 13 MΩ bias resistor and the bias T and are represented by the 170Ω/450 pF RC filter on the biasing line in Supplemental Material Fig. 1. The distributed RC filter has two benefits on the effective electron temperature of our experiment: it provides a high frequency filtering that reduces the polarisation noise as well as thermalisation of the electrons. The copper powder filter is meant to absorb parasitic microwave noise. The line allowing to measure the low frequency response of the junction is filtered by a multipole RC low pass filter, made with a succession of 2 kΩ Nickel-Chromium resistances and 1 nF capacitances to ground. The NiCr resistances were checked in an independent cool-down to change by less than 1%, which allows to calibrate the 13 MΩ resitor in-situ, with a precision better than 1%, which in turn allows us to determine the dc voltage V applied to the tunnel junction. The validity of this callibrattion is confirmed by the quality of the comparison between the observed steps in ∂S I (ν)/∂V and our predictions.

Microwave circuit and calibration

The microwave chain comprises a bias T, two 4-8 GHz cryogenic circulators anchored at the mixing chamber, as well as a 4-8 GHz bandpass filter and a 12 GHz low pass gaussian absorptive filter (see Supplemental Material Fig. 2). These elements are anchored on the mixing chamber and are meant to protect the sample from the back-action noise of the amplifier. The quantitative determination of the detection impedance relies on the detection of the power emitted by the shot noise of the tunnel junction in the high bias regime. We bias the junction at ∼ 1 mV, where DCB corrections are negligible, so that S I = 2eI at frequencies |ν| eV /h 0.5 THz. In order to separate this noise from the noise floor of the cryogenic amplifier, we then apply small variations of the bias voltage and measure the corresponding changes in the measured microwave power with a lock-in amplifier. The conversion of S I into emitted microwave power depends on the environment impedance Z(ν) seen by the tunneling resistance R T . First, only a fraction R 2 T /|R T + Z(ν)| 2 of the current noise is absorbed by the environment. The current noise in the environment has then to be multiplied by Re[Z(ν)] to obtain the microwave power emitted by the electronic shot noise:

S P (ν) = 2eV Re[Z(ν)]G -1 T |Z(ν) + G -1 T | 2 2eV Re[Z(ν)]G -1 T [Re[Z(ν)] + G -1 T ] 2 . ( 5

) The last approximation, Im[Z(ν)]

G -1 T is satisfied with a precision better than 2%. Finally, what is actually detected at room temperature is the amplified microwave power:

S RT P (ν) = 2eV G(ν) Re[Z(ν)]G -1 T [Re[Z(ν)] + G -1 T ] 2 . (6)
Supplemental Material Eq. 6 shows that the extracted Re[Z(ν)] depends on the gain of the microwave chain G(ν), which has to be determined in-situ and independently. To do so, we inserted a 20 dB directional coupler between the sample and the bias Tee, and injected through an independently calibrated injection line, comprising 70 dB attenuation distributed between 4.2 K and the mixing chamber temperature (see Supplemental Material Fig. 2). Both the attenuators and the directional coupler were calibrated at 4.2 K. The resonator is tuned at a resonance well below 4 GHz, so that the entire microwave power is reflected by the sample, allowing to calibrate the gain of the detection chain, albeit doubling the insertion losses of the bias Tee and of the 10 cm microwave cable connecting it to the sample. However, this parasitic contribution can be substracted thanks to independent calibrations. Due to lack of space, we then had to remove one of the circulators for this experiment, which increased the electronic temperature from 16 mK to 25 mK, which explains that we calibrated the detection impedance in an independant run of the experiment.

EXTRACTING THE CURRENT NOISE

We discuss here the possible consequences of the fact that the detection impedance is not negligible compared to the tunneling resistance. More specifically, we show that due to the variations of the tunneling resistance with bias voltage, measuring ∂S P (ν)/∂V is not rigourously equivalent to measuring ∂S I /∂V . However, the error introduced by this approximation can be shown to be negligible.

Due to the non linearity of the tunnel transfer, the power emitted by the junction biased at bias V reads

S P (ν) = Re[Z(ν)] Z T (ν, V ) Z(ν) + Z T (ν, V ) 2 S I (V, ν). (7)
Here Z T (ν, V ) is the differential impedance of the junction, biased at voltage V , at the measurement frequency ν. Supplemental Material Eq. 7 is valid as long as the ac current going through the junction as a consequence of the shot noise is small enough for the response of the junction Z T (ν, V ) to remain in the linear regime. In that case, the modulation of the output voltage of the quadractic detector that we measure is proportional to

∂S P (ν) ∂V = Re[Z(ν)] Z T (ν, V ) Z(ν) + Z T (ν, V ) 2 ∂S I (V, ν) ∂V +S I (V, ν) ∂ ∂V Z T (ν, V ) Z(ν) + Z T (ν, V ) 2 . (8)
We can deduce the expected variations of Z T (ν, V ) with bias voltage from the dc transport properties of the junction via Inserting Eq. 9 in Eq. 8 shows that the associated corrections are negligeable, so that detecting ∂S P (ν) ∂V gives direct access to ∂S I (V,ν) ∂V within a precision better than 1%.

NOTE ON NOISE SYMETRIZATION

Let us add a note about the symetrization of the correlator being measured in our high frequency noise measurements. On one hand the microwave amplifier probes its input voltage [4], which can be written as V in (t) = -Z 0 4π B ω 1/2 (ia in (ω)e -iω(t-x/c) + h.c.)dω, (10) where a in (ω) stands for the destruction operator of a input (right moving) photon at frequency ω, c is the velocity of electromagnetic waves in the transmission line connecting the sample to the amplifier, of characteristic impedance Z 0 = 50Ω, and B the measurement bandwidth. The output voltage reads:

V out (t) = -Z 0 4π B ω 1/2 (ia out (ω)e -iω(t-x/c) + h.c.)dω, (11) with

a out = √ Ga in + √ G -1f † (12)
with G the gain of the amplifier and f † representing the amplifier's noise [4]. Thus the final power measurement ∝< V 2 out > contains, on top of the amplifier's noise, a term proportional to

V 2 in = Z 0 4π B ω a in (ω)a † in (ω) + a † in (ω)a in (ω) dω (13) = Z 0 2π B ω a † in (ω)a in (ω)) + 1 2 dω, (14)
thus containing a power representing the zero point motion of the line, associated to a power ω/2 per unit bandwidth. In other words, the amplifier gives access to the sum of the absorption and emission noise of the input line, unlike quantum detectors [9,10] which allow detecting them separately. However, this does not imply that we measure the electronic current noise of the sample symetrized with respect to frequency: In our experiment, we measure the excess output noise power associated to the dc biasing of our sample. This results in a excess population of the input field < a † in (ω)a in (ω) >, which is itself proportional to the photons leakage rate out of our resonator. Thanks to the circulators anchored at the mixing chamber temperature T hν 0 /k B , the latter is kept very close to its ground state (the maximum average number of photons in the resonator is below 5% for the experiments reported in the paper), and can thus only absorb energy from the electronic current fluctuations (in other words P (hν) P (-hν)). This is why our signal is proportional to the emission noise power S I (ν, V ) and not the noise power symetrized with respect to frequency S sym I (ν, V ) = [S I (-ν, V ) + S I (ν, V )]/2, which also contains the electronic zero point motion.

Figure 1 .

 1 Figure 1. Coulomb Blockade in a normal quantum conductor: (a) A quantum conductor (here a tunnel junction) is voltage biased (V ) through a series impedance Z modeled as a collection of harmonic modes, resulting in inelastic electron tunneling. (b) Experimental set-up: a tunnel junction is dc biased at voltage V , and connected to a SQUID-based resonator presenting a discrete mode at frequency ν0, tuned by varying the magnetic flux φ b threading each SQUID loop. The dc biasing line and 50Ω microwave measurement line are separated by a bias tee, allowing to measure the junction dc differential conductance G(V ) and the emission noise SI (ν0, V ). The measurement line includes an isolator, a cryogenic amplifier with 42 dB gain, a 180 MHz passband filter centered on ν0, and a matched quadratic detector. Temperatures of the different stages are indicated on the right. Bottom illustrations: global view of the sample, with SEM pictures of SQUIDs (top inset) forming the array and of the normal tunnel junction (bottom inset), both from sample 1

Figure 2 .

 2 Figure 2. Characterization of the environment impedance: Designed (left) and measured (right) real part of the impedance Z(ν) of the quarter-wave resonator of sample 1 (see Fig.1) as a function of magnetic flux φ b and frequency ν. The overprinted curve on the right (top scale) shows the resonance at ν0 6GHz for φ b = 0, measured (solid line) and calculated (dotted line). Horizontal dotted white lines indicate -3dB bandwidth used for measuring the shot noise power shown in Fig 3a.

Figure 3 .

 3 Figure 3. Comparison between the measured conductance and noise blockade, and an extension of the dynamical Coulomb blockade theory. Normalized differential conductance G(V ) (top) and current noise spectral density ∂SI(ν0, V )/∂V (bottom). Open circles are experimental data measured at 15 mK. The left panel shows data measured on sample 1, with ν0 = 6 GHz, the center and right panel data measured on sample 2 with ν0 = 6 GHz and 4 GHz, respectively. Solid red lines result from an analytical fit to the data involving series impedance made of three discrete modes shown in insets and the dotted black curve shows the non-interacting, finite frequency shot noise prediction. The green dot-dashed line represents the DCB expression for the current noise density symetrized with respect to frequency.

Z - 1 T

 1 (ν, V ) = e 2hν [I(V + hν/e) -I(Vhν/e)]. (9)

Figure 1 .

 1 Figure 1. More detailled view of the experimental setup. Only circuit components inside the refrigerator are shown in full detail.

Figure 2 .

 2 Figure 2. Characterization of the environment impedance: the microwave chain is calibrated using a vectorial network analyzer, setting the resonance frequency well below 4 GHz.

  Supplemental Material for the article "Dynamical Coulomb Blockade of Shot Noise" Carles Altimiras, * Olivier Parlavecchio, Philippe Joyez, Denis Vion, Patrice Roche, Daniel Esteve, and Fabien Portier

	Service de Physique de l'Etat Condensé (CNRS URA 2464),
	IRAMIS, CEA-Saclay, 91191 Gif-sur-Yvette, France
	(Dated: September 20, 2013)

Blanter20001Y. Blanter and M. Büttiker, Shot noise in mesoscopic conductors," Physics Reports, vol. 336, no. 1-2, pp. 1 -166, 2000. [2] G.-L. Ingold and Y. V. Nazarov, "Charge tunneling rates in ultrasmall junctions," in Single Charge Tunneling (H. Graber and M. H. Devoret, eds.), Plenum Press (New York and London), 1992.


Julien Jacques
email: julien.jacques@univ-lyon2.fr

Christophe Biernacki
email: christophe.biernacki@math.univ-lille1.fr

Model-Based Co-clustering for Ordinal Data

Keywords: coclustering, données ordinales, algorithme SEM-Gibbs model-based coclustering, ordinal data, SEM-Gibbs algorithm

Nous présentons dans ce travail un algorithme de coclustering pour données ordinales. Cet algorithme repose sur le modèle des blocs latents utilisant le modèle BOS (Biernacki et Jacques, 2015, Stat. Comput.) pour données ordinales et un algorithme SEM-Gibbs pour l'inférence. Des expériences sur données simulées illustrent l'efficacité de la méthode d'estimation.

Intoduction

Historically, clustering algorithms are used to explore data and to provide a simplified representation of them with a small number of homogeneous groups of individuals (i.e. clusters). With the big data phenomenon, the number of features becomes itself larger and larger, and traditional clustering methods are no more sufficient to explore such data. Coclustering algorithms have been introduced to provide a solution by gathering into homogeneous groups both the observations and the features. Among numerous algorithms recently developed, model-based approaches [START_REF] Govaert | Co-Clustering[END_REF] have proven their efficiency.

This work focuses on ordinal data, which are categorical data with ordered levels. To the best of our knowledge, the first model-based coclustering algorithm for such data has been proposed by [START_REF] Matechou | Biclustering models for ordinal data[END_REF], relying on the proportional odds model. In this work, we propose a model-based coclustering algorithm relying on the BOS (Binary Ordinal Search) model [START_REF] Biernacki | Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm[END_REF], which has proven its efficiency for modeling ordinal data.

The data set is composed of a matrix of n observations (rows or individuals) of d ordinal variables (columns or features): x = (x ih ) 1≤i≤n,1≤h≤d . For simplicity, the levels of x ih will be numbered {1, ..., m h }, and all m h 's are assumed to be equal:

m h = m (1 ≤ h ≤ d).

Latent block model for ordinal data

Latent block model The latent block model assumes local independence i.e., the n × d random variables x are assumed to be independent once the row partition v = (v ik ) 1≤i≤n,1≤k≤K and the column partition w = (w h ) 1≤h≤d,1≤ ≤L are fixed (note that a standard binary partition is used for v and w):

p(x; θ) = v∈V w∈W p(v; θ)p(w; θ)p(x|v, w; θ) (1)
with (below the straightforward range for i, h, k and are omitted):

• V the set of all possible partitions of rows into K groups, W the set of partitions of the columns into L groups,

• p(v; θ) = ik α v ik k and p(w; θ) = h β w h where α k and β are the row and column mixing proportions, belonging to [0, 1] and summing to 1,

• p(x|v, w; θ) = ihk p(x ih ; µ k , π k ) v ik w h where p(x ih ; µ k , π k ) is the probability of
x ij according to the BOS model [START_REF] Biernacki | Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm[END_REF] parametrized by (π k , µ k ) with the so-called precision parameter π k ∈ [0, 1] and position parameter µ k ∈ {1, ..., m} (detail of p(x ih ; µ k , π k ) is given below),

• θ = (π k , µ k , α k , β ) is the whole mixture parameter.

Ordinal model The BOS model [START_REF] Biernacki | Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm[END_REF], built using the assumption that an ordinal variable is the result of a stochastic binary search algorithm in which e j is the current interval in {1, . . . , n}, and y j the break point in this interval, is defined as follows:

p(x ij ; µ k , π k ) = e m-1 ,...,e 1 m-1 j=1 p(e j+1 |e j ; µ k , π k )p(e 1 ) (2)
where p(e j+1 |e j ; µ k , π k ) = The density (1) can finally be written

p(x; θ) = v∈V w∈W ik α v ik k h β w h ihk p(x ih ; µ k , π k ) v ik w h . ( 3
)
Missing data In the present work, the data x may be also incomplete. We will notice x the set of observed data, x the set of unobserved data and x the set of both observed and unobserved data. The algorithm we propose below will be able to take into account these missing data and to estimate them.

Model inference

The aim is to estimate θ by maximizing the observed log-likelihood

(θ; x) = x ln p(x; θ). (4)
For computational reasons, and EM algorithm is not feasible in that coclustering case (see [START_REF] Govaert | Co-Clustering[END_REF]), thus we opt for one of its stochastic version denoted by SEM-Gibbs [START_REF] Govaert | Estimation d'un modèle à blocs latents par l'algorithme sem[END_REF].

SEM-Gibbs algorithm

The proposed SEM-Gibbs algorithm relies on the EM algorithm used in [START_REF] Biernacki | Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm[END_REF] for the estimation of the BOS model. Starting from an initial value for the parameter θ (0) and for the missing data x(0) , w (0) , the qth iteration of the SEM-Gibbs algorithm alternates the following SE and M steps (q ≥ 0).

SE step

Execute a small number (at least 1) of successive iterations of the following three steps:

1. generate the row partition v

(q+1) ik |x (q) , x, w (q) for all 1 ≤ i ≤ n, 1 ≤ k ≤ K: p(v ik = 1|x (q) , x, w (q) ; θ (q) ) = α (q) k f k (x (q) i. |w (q) ; θ (q) ) k α (q) k f k (x (q) i. |w (q) ; θ (q) )
where

f k (x (q)
i. |w (q) ; θ (q) ) = h p(x

(q) ih ; µ (q) k , π (q)
k ) w (q) h and x (q) i. = {x (q) , x} i denotes observed data {x} i completed by the generation of the unobserved one {x (q) } i at the qth step.

symmetrically, generate the column partition w

(q+1) h |x (q) , x, v (q+1) for all 1 ≤ h ≤ d, 1 ≤ ≤ L: p(w h = 1|x (q) , x, v (q+1) ; θ (q) ) = β (q) g (x (q)
.h |v (q+1) ; θ (q) ) β (q) g (x

(q)

.h |v (q+1) ; θ (q) ) where g (x

(q)

.h |v (q+1) ; θ (q) ) = ik p(x

(q) ih ; µ (q) k , π (q) k ) v (q+1)
ik and x

(q)

.h = {x (q) , x} h

generate the missing data

x(q+1

) ih |x, v (q+1) , w (q+1) following p(x ih |x, v (q+1) , w (q+1) ; θ (q) ) ∝ k α (q) k v (q+1) ik β (q) w (q+1) h k p(x ih ; µ (q) k , π (q) k ) v (q+1) ik w (q+1) h .
M step Estimate θ, conditionally on x(q+1) , v (q+1) , w (q+1) obtained at the SE step (and also conditionally to x), using the EM algorithm of [START_REF] Biernacki | Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm[END_REF].

Choosing the parameter estimation After a burn in period, the final estimation of the discrete parameter µ k is the mode of the sample distribution, and the final estimation of the continuous parameters (π k , α k , β ) is the mean of the sample distribution. It produces a final estimate θ. The missing data estimation x are obtained by the mode of the sample distribution after the burn-in period too.

Choosing the partition The final bi-partition (v, ŵ) is estimated by maximum a posteriori according the parameter estimation.

Likelihood approximation and initialization strategy

Since the observed-data log-likelihood (4) is not tractable, it is approximated by the following harmonic mean [6, chapter 7]:

l( θ; x) ≈ -ln 1 Q -B Q q=B 1 p(x|x (q) , v (q) , w (q) ; θ) (5)
where (x (q) , v (q) , w (q) ) arise independently from p(x (q) , v (q) , w (q) |x; θ) (they are simulated sequentially as in the SE-Gibbs step with Q iterations after a burn in period of length B), and with: p(x|x (q) , v (q) , w (q) ; θ)

∝ ik αv (q) ik k h βw (q) h ihk p(x ih ; μk , πk ) v (q) ik w (q) h (6)
where x ih is x(q) ih if it corresponds to an unobserved data. In order to achieve convergence of the SEM-Gibbs algorithm to a global maximum, multiple initialization is used in practice (typically 10 random initializations), and the result maximizing the approximated likelihood will be selected.

Choice of the number of clusters

In order to select the numbers of clusters, K in rows and L in columns, we propose to adapt to our situation the ICL-BIC criterion developed in [START_REF] Keribin | Estimation and selection for the latent block model on categorical data[END_REF] in the case of coclustering of categorical data:

ICL-BIC(K, L) = log p(x, x, v, ŵ; θ) - K -1 2 log n - L -1 2 log d - KL 2 log(nd) (7)
where x, v, ŵ and θ are the respective estimation of the missing data, row partition, column partition and model parameters obtained at the end of the estimation algorithm.

4 Numerical experiments SEM-Gibbs algorithm validation 30 data sets are simulated according to the following setup: K = L = 3 clusters row and column, d = 12 ordinal variables with m = 5 levels and n = 100 observations. The values of (µ k , π k ) are given by Table 1.

k/ 1 2 3 1 (1,0.9) (2,0.9) (3,0.9) 2 (4,0.9) (5,0.9) (1,0.5) 3 (2,0.5) (3,0.5) (4,0.5) Table 1: Parameter values used for experiments.

We use 50 iterations of the SEM-Gibbs algorithm with a burn-in period of 20 iterations. These numbers seem graphically sufficient to achieve stability of the simulations. Figure 1 illustrates the efficiency of the proposed estimation algorithm, by plotting the coclustering results and the following indicators:

• mu (resp. pi): mean distance between the true µ (resp. π) and its estimated value μ (resp. π): ∆µ

= K k=1 L =1 |µ k -μk |/(KL) (resp. ∆π = K k=1 L =1 |π k - πk |/(KL)),
• alpha (resp. beta): mean distance between the true α (resp. β) and its estimated value α (resp. β): ∆α = K k=1 |α k -αk |/K (resp. ∆β = L =1 |β -β |/L), • ARIr (resp. ARIc): Adjusted Rand Index (ARI) for the row (resp. column) partition.

Ongoing work Further simulation studies (efficiency of the ICL-BIC criterion to select K and L, influence of missing data) and real data analysis will be presented during the conference.

y

  j ∈e j p(e j+1 |e j , y j ; µ, π)p(y j |e j ), p(e j+1 |e j , y j ; µ k , π k ) = π k p(e j+1 |y j , e j , z j = 1; µ k ) + (1 -π k )p(e j+1 |y j , e j , z j = 0), p(e j+1 |y j , e j , z j = 0) = |e j+1 | |e j | I(e j+1 ∈ {e - j , e = j , e + j }), p(e j+1 |y j , e j , z j = 1; µ k ) = I(e j+1 = argmin e∈{e - j ,e = j ,e + j } δ(e, µ k ))I(e j+1 ∈ {e - j , e = j , e + j }), with δ a "distance" between µ and an interval e = {b -, . . . , b + }: δ(e, µ) = min(|µb -|, |µ -b + |) and also p(z j |e j ; π k ) = πI(z j = 1) + (1 -π k )I(z j = 0) and p(y j |e j ) = 1 |e j | I(y j ∈ e j ).

Figure 1 :

 1 Figure 1: An example of data (top left), coclustering results (top right), error in parameter estimation (bottom left) and ARI for the row and column partitions (bottom right).