Yuji Isshiki

Dongzhe Li

Manabu Kiguchi

Tomoaki Nishino

Fabian Pauly
email: fabian.pauly@uni-a.de

Shintaro Fujii
email: fujii.s.af@m.titech.ac.jp

Structural Asymmetry of Metallic Single-Atom Contacts Detected by Current-Voltage Characteristics

Keywords: electrical conductance, metallic atomic contact, break junction, current-voltage characteristics, structural asymmetry

The complex behavior of the simplest atomic-scale conductors indicates that the electrode structure itself is significant in the design of future nanoscale devices. In this study, the structural asymmetry of metallic atomic contacts formed between two macroscopic Au electrodes at room temperature was investigated. Characteristic signatures of the structural asymmetries were detected by fast current-voltage (I-V) measurements with a time resolution of approximately 100 µs. Statistical analysis of more than 300,000 I-V curves obtained from more than 1,000 contactstretching processes demonstrates that the current rectification properties are correlated with the conductance of the nanocontacts. A substantial suppression of the variation in current rectification was observed for the atomic contacts with integer multiples of the conductance quantum.

Statistical analysis of the time-resolved I-V curves revealed that the current rectification variations increased significantly from 500 µs onwards before the breakage of the atomic contacts. Atomistic ab initio simulations of the stretching processes and corresponding I-V characteristics confirmed the magnitude of the rectification and related it to the structural asymmetries in the breakdown process of the junctions. This study provides a better understanding of charge transport properties in extremely sensitive quantum point contacts, where the geometric and the electronic structure at the atomic scale are closely related each other.

INTRODUCTION

Nanocontacts fabricated using break junction and scanning probe methods have helped scientists investigate concepts in charge transport at the atomic scale [START_REF] Agraït | Quantum Properties of Atomic-Sized Conductors[END_REF][START_REF] Van Ruitenbeek | Conductance Quantization in Metallic Point Contacts[END_REF]. In the last few decades, experimental studies on charge transport across metallic atomic contacts have revealed unique quantum properties, such as electronic shell effects [START_REF] Yanson | Observation of Shell Structure in Sodium Nanowires[END_REF][START_REF] Mares | Observation of Electronic and Atomic Shell Effects in Gold Nanowires[END_REF], parity oscillations [START_REF] Smit | Observation of a Parity Oscillation in the Conductance of Atomic Wires[END_REF], conductance fluctuations [START_REF] Ludoph | Evidence for Saturation of Channel Transmission from Conductance Fluctuations in Atomic-Size Point Contacts[END_REF], and quantum interference [START_REF] Ludoph | Evidence for Saturation of Channel Transmission from Conductance Fluctuations in Atomic-Size Point Contacts[END_REF]. Most of the experiments were conducted under cryogenic conditions. At cryogenic temperatures and fixed electrode separation, an atomic contact is trapped in a specific geometrical configuration. In contrast, large variations in contact geometry can arise from thermal fluctuations at room temperature [START_REF] Muller | Quantization Effects in the Conductance of Metallic Contacts at Room Temperature[END_REF][START_REF] Rubio | Atomic-Sized Metallic Contacts: Mechanical Properties and Electronic Transport[END_REF][START_REF] Kizuka | Atomic Process of Point Contact in Gold Studied by Time-Resolved High-Resolution Transmission Electron Microscopy[END_REF][START_REF] Kizuka | Atomic Configuration and Mechanical and Electrical Properties of Stable Gold Wires of Single-Atom Width[END_REF][START_REF] Dreher | Structure and Conductance Histogram of Atomic-Sized Au Contacts[END_REF]. In this context, the linearity and nonlinearity of current versus bias voltage (I-V) characteristics of metallic atomic contacts have been investigated under thermal fluctuation conditions at room temperature [START_REF] Yoshida | Nonlinear Current-Voltage Curves of Gold Quantum Point Contacts[END_REF][START_REF] Costa-Krämer | Conductance Quantization in Nanowires Formed between Micro and Macroscopic Metallic Electrodes[END_REF]. Hansen et al. developed a method for fast charge transport measurements to capture the transport properties of metastable states of single-atom contacts at room temperature [START_REF] Hansen | Fast and Accurate Current-Voltage Curves of Metallic Quantum Point Contacts[END_REF]. The method detects the I-V characteristics, and has been successfully applied to reveal the linear I-V behavior of Au contacts at room temperature [START_REF] Hansen | Fast and Accurate Current-Voltage Curves of Metallic Quantum Point Contacts[END_REF][START_REF] Hansen | Current-Voltage Curves of Gold Quantum Point Contacts Revisited[END_REF]. Regarding the linearity of the I-V behavior, Tewari et al. reported shot noise measurements at liquid He temperature [START_REF] Tewari | Anomalous Nonlinear Shot Noise at High Voltage Bias[END_REF]. They revealed that the differential conductance (i.e., dI/dV) of Au atomic contacts can be both antisymmetric and asymmetric. The asymmetric character was attributed to electron backscattering from defect sites close to the atomic contact [START_REF] Tewari | Anomalous Nonlinear Shot Noise at High Voltage Bias[END_REF]. Therefore, the I-V characteristics of metallic atomic contacts are sensitive to their structural symmetry. Thus, this study aims to identify the effects of the structural symmetry of metallic atomic contacts at room temperature, which has not yet been achieved using conventional transport measurements. To this end, a fast I-V measurement setup with a high temporal resolution of approximately 100 µs was developed using the field-programmable gate array (FPGA) technology to sample a number of metastable states of Au atomic contacts at room temperature.

We demonstrate that Au atomic contacts with a conductance equal to integer multiples of the conductance quantum, G0 = 2e 2 /h, exhibit a suppressed asymmetry in the I-V characteristics, and hence demonstrate reduced rectification effects in comparison to that of the nanocontacts with a noninteger conductance. Time course analysis of the Au atomic contacts at breakdown reveals that the structural asymmetry considerably increases starting at 500 µs before the breakage of the contacts. Thus, the Au atomic contacts were subjected to structural fluctuations, in which the contacts became asymmetric. By combining experiments and charge transport simulations, the rectifying properties of nanojunctions are related to their structural asymmetry. Consequently, this study reveals the dynamic behavior of Au atomic contacts formed during the stretching process at room temperature at a temporal resolution of approximately 100 µs.

RESULTS AND DISCUSSION

Conductance and Current Rectification Ratio of Au Atomic Contacts

Au atomic contacts were created by repeatedly forming and breaking the contact between the Au(111) substrate and Au STM tip under a bias voltage of 100 mV (Figure 1a; see also Supporting Information 1 and 2). The electrical current was recorded while forming and breaking the nanocontact. When the electrical conductance of the contact at a voltage of 100 mV (G = I/V) decreased to below 4G0 in the breaking process, the tip position was fixed, and the bias voltage was cycled between -92 and +92 mV (Figures 1b and1c). One cycle of the voltage sweep was set to 100 μs. After the I-V measurements in the self-breaking process, a nanocontact with a conductance above 4G0 was immediately reformed by pushing the Au STM tip toward the Au(111) substrate. Measurements in the self-breaking process mode of the atomic contacts were repeated until a statistically significant dataset was obtained; a total of 309,600 I-V curves were collected, corresponding to 2,464 self-breaking events. Figures 1d and1e present the time evolution of the conductance and rectification behavior of an Au contact during a self-breaking process. The conductance in the positive and negative bias regimes (i.e., Gp and Gn) was calculated by linearly fitting each I-V curve in the bias range from -70 to 0 mV and from 0 to +70 mV, respectively.

The rectification behavior was analyzed by calculating the current rectification ratio (RR), which is defined as RR = Gp/Gn. Here, Gp and Gn are the calculated conductances at the positive and negative bias voltages, respectively. In the time evolution of the self-breaking process, the conductance features plateaus around 2.5G0, 2.0G0, 1.5G0, and 1.0G0 at both positive and negative bias voltages, which are interconnected by sudden conductance jumps (Figure 1d). Based on our measurement resolution, the conductance jumps occur within time intervals of 0.10 ms. These results indicate a diameter of only a few atoms at the narrowest cross section of the contact and spontaneous structural rearrangements during the self-breaking process [START_REF] Yanson | Formation and Manipulation of a Metallic Wire of Single Gold Atoms[END_REF][START_REF] Tsutsui | Formation and Self-Breaking Mechanism of Stable Atom-Sized Junctions[END_REF][START_REF] Tsutsui | Heat Dissipation in Quasi-Ballistic Single-Atom Contacts at Room Temperature[END_REF]. At the end of the self-breaking process, the conductance drops from approximately 1.0G0 to the noise floor within 9.5 ms, which signals the breakdown of the Au single-atomic contact. Similar self-breaking behavior has been reported in the studies by Muller et al. [START_REF] Muller | Quantization Effects in the Conductance of Metallic Contacts at Room Temperature[END_REF] and Hansen et al. [START_REF] Hansen | Fast and Accurate Current-Voltage Curves of Metallic Quantum Point Contacts[END_REF], in which Au contacts with conductances below 10G0 spontaneously broke at room temperature. While the Au atomic contacts can be held stable for rather long times on the order of an hour at low temperature, Au contacts spontaneously break at room temperature on a time scale of milliseconds owing to the diffusion of the atoms making up the contact [START_REF] Agraït | Quantum Properties of Atomic-Sized Conductors[END_REF]. The corresponding diffusion of Au atoms and the deformation process of Au atomic contacts at room temperature have been observed in nanoporous Au by time-resolved high-resolution transmission electron microscopy [START_REF] Liu | Time-Resolved Atomic-Scale Observations of Deformation and Fracture of Nanoporous Gold under Tension[END_REF]. Timeresolved high-resolution transmission electron microscopy imaging on the time scale of milliseconds (i.e., a rate of 40 frames per second) revealed that a combination of dislocation-related plastic deformations and slip-activated surface diffusion of Au atoms with a velocity of ca. 30 nm/s leads to stepwise atomistic structural changes, and finally the failure of an Au contact with an initial diameter of approximately 5 nm. The changes in conductance and RR of the Au atomic contacts, measured in the sub-millisecond time domain (Figures 1d and1e), are attributed to structural changes in the Au contacts caused by the deformation and surface diffusion. Because the time scale of the surface diffusion of Au atoms is faster than the time resolution of this study, we expect that the measurements represent an average over different geometric configurations. By not considering the I-V characteristics of Au atomic contacts with kinks, we explored Au atomic contacts before and after major structural rearrangements occurred, as discussed below thoroughly.

Figure 1e shows the time evolution of log10(RR), where fine variations in the current rectification behavior are apparent. In addition, the sign of log10(RR) is an indicator of the direction of asymmetry (positive or negative asymmetry), and the absolute value of log10(RR) characterizes the degree of asymmetry. Figures 2a-d present further examples of conductance plateaus with a G of approximately 1.0G0 and 1.5G0 [START_REF] Muller | Quantization Effects in the Conductance of Metallic Contacts at Room Temperature[END_REF][START_REF] Dreher | Structure and Conductance Histogram of Atomic-Sized Au Contacts[END_REF][START_REF] Mosso | Heat Transport through Atomic Contacts[END_REF], and their corresponding log10(RR). The conductance plateaus with integer multiples of the conductance quantum exhibit only small variations in G (Figure 2a) and log10(RR) (Figure 2b). In sharp contrast, marked fluctuations of G and log10(RR) are visible during the conductance plateaus with half-integer multiples of G0 (Figures 2c and2d). This result agrees with the observations of quantum oscillations in the thermopower, which show reduced thermopower values at integer multiples of the conductance quantum [START_REF] Evangeli | Quantum Thermopower of Metallic Atomic-Size Contacts at Room Temperature[END_REF]. This can be explained by a particularly robust electron transmission formed by fully open transmission eigenchannels [START_REF] Mosso | Heat Transport through Atomic Contacts[END_REF][START_REF] Pauly | Molecular Dynamics Study of the Thermopower of Ag, Au, and Pt Nanocontacts[END_REF].

Statistics of Conductance and Current Rectification Ratio

We collected 309,600 I-V curves to investigate the statistics of the conductance and current RR.

Moreover, 31% of the 309,600 I-V curves displayed kinked shapes. Such I-V curves were automatically removed from the data set (Supporting Information 3) used for the statistical analysis (Figures 3 and4), because the kinked behavior must be owing to a major dynamic structural instability of the Au contacts during a single I-V measurement cycle. The kinked behavior tends to appear at higher bias conditions within the I-V measurement, and after the observation of a kink, the conductance changed by approximately a conductance quantum (Supporting Information 3).

Figure 3a shows a two-dimensional (2D) histogram of log10(RR) versus G, which was constructed using 212,279 I-V curves without kinks (i.e., 69% of 309,600). It reveals three regions with high intensity counts located near 1G0, 2G0, and 3G0, in which log10(RR) is symmetrically distributed near 0. This reflects an energetic preference for Au contact geometries with a quantized conductance, and a symmetric contact structure formed on average between the Au tip and substrate. The highest-intensity regions with quantized conductance values of 2G0 and 3G0 are slightly shifted toward lower conductance values of 1.9G0 and 2.7G0, respectively. This conductance shift agrees with the results obtained in a previous study [START_REF] Agraït | Quantum Properties of Atomic-Sized Conductors[END_REF], in which the origin of the shift was attributed to electron backscattering on defects near the contact. In addition, we investigated the temperature dependence of the conductance and current RR statistics (Supporting Information 4). Heating to approximately 320 K suppressed the formation of geometries, exhibiting conductances of more than 1G0 (i.e., 2G0 and 3G0), which confirms previously reported transport measurement results at a fixed bias voltage [START_REF] Mosso | Heat Transport through Atomic Contacts[END_REF][START_REF] Cui | Quantized Thermal Transport in Single-Atom Junctions[END_REF][START_REF] Mosso | Quantized thermal conductance in metallic heterojunctions[END_REF].

To analyze the relationship between the conductance and structural symmetry of the metallic nanocontacts, log10(RR) histograms were constructed using a selected dataset with conductance values in the ranges of (1.0 ± 0.036)G0, (1.5 ± 0.036)G0, and (1.9 ± 0.036)G0 (Figure 3b). The corresponding log10(RR) histograms were fitted using the Lorentzian function. It can be observed that the full width at half maximum (FWHM) of the Lorentzian fit w is suppressed when the conductance values are multiples of the conductance quantum. This conductance-dependent suppression of w is more clearly visible in the conductance versus w plot (Figure 3c), in which w oscillates with the conductance, and the minima of w are clearly located at multiples of G0.

Dynamics of the Rectification Ratio before Rupture of Au Atomic Contacts

Figures 4a and4b show an example of the last conductance plateau and corresponding log10(RR) before the breakage of an Au atomic contact with a conductance of 1G0. Immediately before the breakdown, the conductance gradually changes by approximately 30%, and finally decreases rapidly to the noise floor within 100 µs. In contrast to the conductance, log10(RR) shows remarkable variations with time. Initially, they are near 0; however, the magnitude of log10(RR) drastically increases shortly before rupture. To obtain statistical information on the dynamic behavior of the conductance and log10(RR), data analysis was performed as follows: First, the contact-breakdown time (t = t0) was defined as the time when the conductance dropped by more than 0.5G0 to the noise floor in the dataset of the nonkinked 212,279 I-V curves that correspond to 1,240 self-breaking events. Because the time resolution of conductance is 0.05 ms and that of the RR is 0.1 ms, there may be no data for the RR at t = t0. For example, the conductance shows 1G0 at t -t0 = 0 and decreases to 0.3G0 or less at t -t0 = 0.05 ms, as shown in Figure 4a, indicating that the Au atomic contact breaks at t = t0. Subsequently, the RR immediately before the breakdown is defined at t -t0 = -0.05, -0.15, -0.25 ms etc., as shown in Figure 4b. Among the 1,240 selfbreaking events, 1,088 realizations exhibited large conductance drops with ∆G>0.5G0.

Subsequently, log10(RR) histograms were generated for the dataset of 1,088 self-breaking events at t -t0 = -0.05, -0.15, and -0.25 ms (Figure 4c). The log10(RR) histograms illustrate a symmetric peak around 0. In Figure 4d, w of the symmetric peak is plotted as a function of time (i.e., t-t0).

The FWHM w starts to increase near -0.5 ms, finally growing to more than 1.5 times the value at -0.05 ms. Thus, statistical analysis of the Au atomic contact-breakdown dynamics indicates that the growing variation in log10(RR) signals the electrical breakdown already 0.5 ms in advance.

Analysis of the evolution of log10(RR) suggests that 18% of the self-breaking processes (i.e., 453/2,464) exhibit a directional evolution of the asymmetry, as indicated by a monotonically increasing (decreasing) value of log10(RR). In 61% of the breaking processes (i.e., 1,499/2,464), log10(RR) instead behaves randomly and jumps between positive and negative values (Supporting Information 5).

Origin of the Rectification Behavior

Theoretical simulations based on density functional theory [START_REF] Giannozzi | QUANTUM ESPRESSO: a Modular and Open-Source Software Project for Quantum Simulations of Materials[END_REF] were performed using three model structures with increasing degrees of asymmetry to explore the origin of the observed conductance and rectification behaviors of the Au atomic contacts, and to relate them to junction geometries, as shown in the top panels of Figure 5 (see also Supporting Information 6). Models 1, 2, and 3 were constructed with 4×4|3×3|3×3|4×4, 4×4|2×2|3×3|4×4, and 4×4|1×1|2×2|3×3|4×4 layers as a dynamic region, which was energetically optimized during the stretching process until the junction broke (Figure 5). Subsequently, we computed the electronic structure in the scattering region selfconsistently for different bias voltages and determined the energy-and bias-dependent transmission function, T(E,V), using the nonequilibrium Green's function formalism [START_REF] Brandbyge | Density-Functional Method for Nonequilibrium Electron Transport[END_REF][START_REF] Papior | Improvements on Non-Equilibrium and Transport Green Function Techniques: The Next-Generation Transiesta[END_REF],

where the scattering region is embedded into the left and right semi-infinite electrodes via self-energies. We employed T(E,V) to obtain the charge current using the Landauer-Büttiker formula as follows [START_REF] Brandbyge | Density-Functional Method for Nonequilibrium Electron Transport[END_REF][START_REF] Papior | Improvements on Non-Equilibrium and Transport Green Function Techniques: The Next-Generation Transiesta[END_REF]:

𝐼𝐼(𝑉𝑉) = 2𝑒𝑒 ℎ � 𝑑𝑑𝑑𝑑 +∞ -∞ 𝑇𝑇(𝑑𝑑, 𝑉𝑉)[𝑓𝑓(𝑑𝑑, 𝜇𝜇 L ) -𝑓𝑓(𝑑𝑑, 𝜇𝜇 R )], (1)
where f(E, µX) is the Fermi-Dirac distribution function and µX is the electrochemical potential of the electrode X = L, R. Finally, RR = I(+100 meV)/I(-100 meV) was evaluated by calculating the nonequilibrium charge currents at the positive and negative bias voltages of V = ±100 meV. The effect of the applied voltage V was included by shifting the electrochemical potential of the left (substrate) and right (STM tip) electrodes by +V/2 and -V/2, respectively. Let us note that because of the self-consistent electronic structure calculation, I(V) is insensitive to how the bias voltage is applied, and we might also have chosen to modify, for instance, only the electrochemical potential of the substrate or, alternatively, those of the tip.

Figure 5 shows the calculated zero-bias conductance and nonequilibrium RR during the contact-stretching process. In contrast to the experiment, the stretching curves are shown as a function of pulling distance, not time. Because the mechanical properties of the flexible electrodes are expected to deviate from those of the experiment owing to the small number of dynamic atoms in the simulations, we considered the theoretically accessible electrode pulling distance.

Remarkably, variations in the conductance and log10(RR) of the model structures exhibit a reasonable quantitative agreement with the experimental results, as shown in Figures 1e and2c. In particular, a relatively pronounced rectification effect is observed for model 3 at the breaking point, where the junction geometries become highly asymmetric. The agreement between the experimental and theoretical results suggests that rectification is associated with the asymmetric junction structure. In addition, the comparison between the experimental and theoretical results may be facilitated by the structural simplicity of the atomic contact at the breaking point because only a single atom or a short chain of atoms is expected to bridge the gap between the two Au electrodes [START_REF] Dreher | Structure and Conductance Histogram of Atomic-Sized Au Contacts[END_REF].

Figure 6 illustrates the energy-dependent transmission functions at ±100 mV for configurations close to the breaking point to achieve a profound insight into the rectification behavior (Supporting Information 7). Depending on the degree of junction asymmetry, the simulated transmissions deviate with respect to the polarity of the applied bias voltage, which results in the rectification property. A slight difference is found in the transmission functions for the positive and negative bias voltages of the rather symmetric models 1 and 2, while a significant difference in T(E,±100 mV) is observed for the asymmetric junction model 3. Therefore, the RR for model 3 is approximately one order of magnitude larger than those of models 1 and 2.

Junctions with a quantized zero-bias conductance exhibit fully transmitting eigenchannels at the Fermi energy. Thus, the transmissions of eigenchannels reach their maximum value of 1, and the full transmission T(E, V = 0) consequently also exhibits a maximum as the sum of the eigenchannel transmissions. Such an extremal curve is highly robust to any changes [START_REF] Evangeli | Quantum Thermopower of Metallic Atomic-Size Contacts at Room Temperature[END_REF][START_REF] Pauly | Molecular Dynamics Study of the Thermopower of Ag, Au, and Pt Nanocontacts[END_REF],

including the application of a relatively small bias voltage of |V|≤100 mV. Therefore, variations in the RRs are expected to be small at zero-bias conductance values close to the quantized values nG0, while they will grow for junctions with intermediate conductances.

CONCLUSIONS

We demonstrated that the I-V characteristic of the Au atomic contacts in the low-bias regime (below 0.1 V) depends on the symmetry of the geometric structure. The I-V curves of contacts at integer multiples of the conductance quantum exhibited an enhanced symmetry compared to the curves at noninteger conductance. Statistical analysis of the I-V characteristics revealed that the variation in the current rectification was suppressed in the contacts with quantized conductances.

This variation demonstrated quantum oscillations as a function of conductance, which are reminiscent of thermopower and conductance fluctuations [START_REF] Evangeli | Quantum Thermopower of Metallic Atomic-Size Contacts at Room Temperature[END_REF][START_REF] Ludoph | Conductance Fluctuations as a Tool for Investigating the Quantum Modes in Atomic-Size Metallic Contacts[END_REF]. Further investigations on the time-dependent I-V behavior revealed that the variation of the current rectification increased significantly during the self-breaking process starting from approximately 0.5 ms before contact rupture. Therefore, the Au atomic contacts are subject to structural fluctuations before the electrical instability point, which renders the contacts asymmetric. Consequently, this study provides a profound insight into the relationship between the geometric and electronic structures of metallic atomic contacts through a fast current-voltage measurement setup under room temperature conditions.

EXPERIMENTAL METHODS

Experiment. An Au(111) surface was prepared by thermal deposition of Au on mica at elevated temperatures under high-vacuum conditions. The employed Au(111) was subsequently flameannealed and washed with ethanol (Kanoto Chemical, purity > 99.5%) before use. A freshly cut Au wire (Nilaco, diameter ~0.3 mm, purity>99%) was used as the tip. A commercially available scanning tunneling microscopy (STM) apparatus was used with a signal access module III (Bruker, Multimode STM), in which the current amplifier and piezo driver were replaced by an SR570 (Stanford Research System) and E-663 (Physik Instrumente), respectively. The feedback loop for the STM operation was controlled using a custom-made program based on an FPGA (National Instrument, PCIe 7852R) with a sampling rate of up to 1 MHz (Supporting Information 1). Gold atomic contacts were created by repeatedly forming and breaking the contact between the Au(111) substrate and Au STM tip with a tip displacement velocity of 7.5 nm/s under a bias voltage of 100 mV (Figure 1a). The electrical current was recorded while forming and breaking the nanocontact.

When the electrical conductance of the contact at 100 mV (G = I/V) decreased to below 4G0 in the breaking process, the tip position was fixed and the bias voltage was swept between -92 and +92

mV. One cycle of the voltage sweep was set to 100 μs. The voltage cycling scheme was continued n × 100 times, where n is a positive integer until the atomic contact is broken. In addition, the conductance decreases to below 0.7G0 within the nth voltage cycling process (Figures 1b and1c; see also Supporting Information 1). After the self-breaking process of the atomic contact at the constant tip-substrate separation, a nanocontact with a conductance above 4G0 was immediately reformed by pushing the Au STM tip toward the Au(111) substrate. Measurements of the selfbreaking processes of the atomic contacts were repeated until a statistically significant dataset was obtained. We collected 309,600 I-V curves, which corresponded to 2,464 self-breaking events.

Calculations. Atomistic ab initio simulations of Au nanowires were carried out within the density functional theory (DFT) framework, as implemented in the plane-wave software package QUANTUM ESPRESSO [START_REF] Giannozzi | QUANTUM ESPRESSO: a Modular and Open-Source Software Project for Quantum Simulations of Materials[END_REF]. The Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional [START_REF] Perdew | Generalized Gradient Approximation Made Simple[END_REF] and ultra-soft pseudopotentials [START_REF] Vanderbilt | Soft Self-Consistent Pseudopotentials in a Generalized Eigenvalue Formalism[END_REF] were used. Energy cutoffs of 30 and 300 Ry were employed for the wave functions and charge density, respectively. Supercell Au structures were obtained by extending the Au(111) unit cell in three dimensions with a 4 × 4 periodicity in the xy plane. The Au contacts were allowed to relax until the atomic forces fell below 10 -3 Ry/Bohr.

The Γ-point approximation was used for geometry optimization. Subsequently, based on the determined geometries, ab initio charge transport properties were computed using the software package TRANSIESTA [START_REF] Brandbyge | Density-Functional Method for Nonequilibrium Electron Transport[END_REF][START_REF] Papior | Improvements on Non-Equilibrium and Transport Green Function Techniques: The Next-Generation Transiesta[END_REF] and its "post-processing" code TBTRANS and SISL [32], which employs the nonequilibrium Green's function formalism combined with DFT. In these separate DFT calculations for transport, Troullier-Martins norm-conserving pseudopotentials, the PBE exchange-correlation functional, the double-zeta polarization basis, and an energy cutoff for the real-space mesh of 250 Ry were used. A 10 × 10 × 1 k-point mesh was found to be sufficient to obtain well-converged transmission functions. We resolved the transmission functions on an energy mesh with a resolution of 0.02 eV.

ASSOCIATED CONTENT Supporting Information

The Supporting Information is available free of charge at

The Supporting Information discusses the temporal resolution of the current measurement setup, examples of I-V measurements with time evolution, analysis of the I-V curves of Au atomic contacts, temperature dependence of conductance and RR, evolution of junction asymmetry, distribution of elongation lengths until rupture, and DFT simulations of Au atomic contacts.

FiguresFigure 1 .

 1 Figures

  Figure 1. (a) Schematic of the experimental setup. (b) Applied bias voltage (V) and example of

Figure 2 .

 2 Figure 2. Examples of the time evolution of (a,c) conductance and (b,d) rectification ratio (RR).

Figure 3 .

 3 Figure 3. (a) Contour plot of counts for corresponding conductance and logarithmically scaled RR

Figure 4 .

 4 Figure 4. Example of the time evolution of (a) conductance and (b) RR before breakage of an Au

Figure 5 .

 5 Figure 5. Evolution of the zero-bias conductance and nonequilibrium RR on a logarithmic

Figure 6 .

 6 Figure 6. Energy-dependent transmission functions at -100 (black) and +100 meV (red) for (a)

ACKNOWLEDGMENTS

This work was financially supported by Grants-in-Aid for Scientific Research (No. 20K05245), JST CREST (No. JP-MJCR18I4), JSPS Fellows (19J21186) from MEXT Japan, and the Precise Measurement Technology Promotion Foundation (PMTPF). This work was performed using HPC resources from CALMIP (Grant 2021-P21008).

All authors contributed in writing the manuscript. All authors approved the final version of the manuscript.

Notes

The authors declare no competing financial interest.

AUTHOR INFORMATION Author Contributions


Tereza Nečasová

Ninon Burgos

David Svoboda
email: svoboda@fi.muni.cz

Validation and evaluation metrics for medical and biomedical image synthesis

Keywords: Image synthesis, Image simulation, Similarity, Distance Metrics, Image Datasets

published or not. The documents may come

Validation and evaluation metrics for medical and biomedical image synthesis

Introduction

Regardless of the application domain, each newly introduced method requires a proper validation procedure. When designing a segmentation algorithm, for example, one is excepted to provide some tests showing how the algorithm performs and whether the results are better compared with state of the art methods. In the field of image synthesis, one can legitimately ask about the explained variability in the synthetic data, whether the generated data look realistic, and whether they are sufficiently similar to their real counterpart. However, no standard pipeline currently exists for synthetic image quality assessment (see Fig. 1). Some authors evaluate their computer-generated data with the statement of having visually similar results, some attempt to measure the similarity using expert's assessment or various techniques. Without any broader comparison of synthesized data, there is no warranty regarding their quality. The quality assessment of the synthesized images should be an integral part of the image synthesis process. The user of the proposed synthetic image data should be assured about their plausibility before using them in further steps.

This chapter includes an overview of the validation methods applied in the field of biomedical and medical image synthesis. Three strategies can be considered, and possibly combined, when evaluating the quality of generators producing synthetic images. In the first case, the data are validated using the strength of expert knowledge to support the plausibility of the generated images. The second class is focused on synthetic images paired with their real counterparts. The third class of evaluation is done by comparing the characteristics of the whole dataset of real images against the characteristics in the whole dataset of synthetic images. Although the assessment in the first class can be further evaluated, it is rather more qualitative than quantitative compared with the two other classes. The main approaches of each strategy will be described in the following subsections.

Expert Knowledge

When developing a new framework responsible for image synthesis where the outputs are expected to resemble the real images, the most expected approach is to report the visual inspection of the synthetic data.

Rating based on the visual plausibility of synthetic data

Figure 1: Processes connected with the validation of synthetic data. Synthetic images obtained from the simulation framework should be compared with real images with a clear conclusion given by a validation method (see in the yellow box) whether they are sufficiently similar or not.

images depicting Pap Smear specimens [START_REF] Malm | Simulation of bright-field microscopy images depicting Pap-smear specimen[END_REF]. The sensitivity relates to the ability of accurate detection of real images as being real and the specificity to the ability of accurate detection of synthetic images as synthetic. These measures are computed from the number of true positives/real images (TP), true negatives/synthetic images (TN), false positives (FP) and false negative cases (FN) as: sensitivity = T P T P + F N

(1) specificity = T N F P + T N .

Both values range from 0 to 1 (=ideal ). In this work, six experts from different research fields were asked to recognize the origin of a set of real and synthetic images. The experts were showed a tightly cropped, randomly selected, view of the real and synthetic scenes projected through lightly frosted glass, to account for limitations of existing display devices. They were showed only for two seconds to imitate the real conditions for assessment. The results of the experts' classification is reviewed and compared to a random classification in Table 1.

A similar validation approach can be found in [START_REF] Calimeri | Biomedical data augmentation using generative adversarial neural networks[END_REF], where the experts were also asked to recognize real or synthetic images resulting from a generator of magnetic resonance imaging (MRI) brain slices. They also compute TP, TN, FP, FN to derive accuracy, precision, sensitivity and F1-score.

In [START_REF] Apou | Synthesizing whole slide images[END_REF], the synthetic whole slide histological images accompanied by the reference images are evaluated by observing the image data at different levels of magnification. The goal is to discover in which level of detail (magnification) the synthetic image is not looking realistic anymore. The measurement of realism of synthetic images is evaluated as follows:

1. open a synthetic image in a software slide viewer on a standard computer screen;

2. decrease the magnification until it "feels realistic";

3. slowly increase the magnification until it "feels wrong"; 4. the current magnification is then written down as the highest realistic magnification for the current image.

The highest realistic magnification is regarded as the score that evaluates the quality of the image simulation. In order the keep some biological correctness, only experts were asked to assess the data plausibility.

The five-point Likert scale (very poor / poor / satisfactory / good / very good ) has been used to assess the quality of synthetic tissue microscopy images generated from pre-segmented Haematoxylin and Eosin images of brain tumours [START_REF] Glotsos | Image quilting and wavelet fusion for creation of synthetic microscopy nuclei images[END_REF]. Here, the expert histopathologists were asked to classify the images according to the following three criteria:

1. similarity to real-world tissue microscopy images, 2. reproducibility of nuclei morphometry and 3. reproducibility in nuclei texture.

However, the particular results of the expert have not been reported.

The Likert scale has also been applied in [START_REF] Vitale | Improving realism in patientspecific abdominal ultrasound simulation using cyclegans[END_REF], where 21 experts in ultrasound were asked to rank (fake, rather fake, cannot decide, rather real, real) the Table 2: Results of the experts' evaluation in [START_REF] Kovacheva | A model of the spatial tumour heterogeneity in colorectal adenocarcinoma tissue[END_REF]. The values report the average evaluation of the appearance of synthetic images by three pathologists. Healthy (H), well differentiated (WD), moderately differentiated (MD), and poorly differentiated (PD) images were evaluated at magnifications 20× and 40×. (1 = Not realistic at all, 5 = Very realistic, '-' means feature is not relevant)

realism of the simulated ultrasound scans. The scans were randomly presented to the experts without giving them any hint whether observing the real or simulated one. The time spent during the analysis of each particular image was also taken into account when assessing the experts' reliability. The objective of the work of Gong et al. [START_REF] Gong | Deep learning enables reduced gadolinium dose for contrast-enhanced brain MRI[END_REF] was to reduce gadolinium dose in contrast-enhanced brain MRI. Two neuro-radiologists were asked to assess the quality of the post-contrast MR images (low-dose, synthesized fulldose and true full-dose). They rated the general image quality, suppression of aliasing/motion artifacts and degree of enhancement compared against precontrast MR images using a five-point Likert scale ranging from 1 (poor) to 5 (excellent).

In [START_REF] Kovacheva | A model of the spatial tumour heterogeneity in colorectal adenocarcinoma tissue[END_REF], a model of healthy and cancerous colonic crypt micro-environment was proposed and successfully implemented to show the ability to control cancer grade, cellularity, cell overlap ratio, and image resolution. Here, the histopathologists were asked to grade the quality of synthesis on a scale from 1 to 5 (5 = very realistic). The averages of the grades for a combination of four types of cell objects, architecture and number of crypts, two different magnifications, and four differentiation levels are reported in Table 2.

When comparing the real and synthetic (also called as fake) images, one can also use a Visual Turing Test (VTT) [START_REF] Geman | Visual turing test for computer vision systems[END_REF]. This test is a procedure during which a stochastic sequence of binary questions is generated and given to some respondents. When applied in the field of medical image analysis, the experts are asked to distinguish between the real and synthetic images in a sequence. VTT was applied for example in [START_REF] Han | Gan-based synthetic brain mr image generation[END_REF] for the validation of generated brain MR images. In another work [START_REF] Han | Synthesizing diverse lung nodules wherever massively: 3d multi-conditional gan-based ct image augmentation for object detection[END_REF], the authors tested the quality of synthesized lung nodules for X-ray computed tomography (CT) image augmentation potentially used for object detection. Chuquicusma et al. [START_REF] Chuquicusma | How to fool radiologists with generative adversarial networks? a visual turing test for lung cancer diagnosis[END_REF] performed VTT on images of malignant and benign lung nodules for a computer-aided diagnosis system generated by a deep convolutional generative adversarial network (DC-GAN). In the context of anomaly detection, Schlegl et al. [START_REF] Schlegl | F-AnoGAN: Fast unsupervised anomaly detection with generative adversarial networks[END_REF] used the VTT to quantify the quality of the generated normal images.

Rating based on the usability of synthetic data

Experts have been asked to not only rate the plausibility of synthetic data but also their usability.

A simulation framework called SIMCEP [START_REF] Lehmussola | Computational framework for simulating fluorescence microscope images with cell populations[END_REF] forms a cornerstone in image synthesis dedicated to fluorescence microscopy. To validate the generated images, mediated experts were asked to use four different image processing tools developed for automated image cytometry (specifically for cell enumeration). Five sets of ten images, each containing 1000 cells and different levels of overlap were analyzed with each out of four tools developed by independent research groups and the results were compared in a plot. The tools gave similar results supporting the expectation that for worse conditions such as overlapping cells, the number of enumerated cells will be lower.

In order to perform a clinical assessment of a method using a generative adversarial network (GAN) to synthesize standard-dose PET images from lowdose ones [START_REF] Ouyang | Ultra-low-dose pet reconstruction using generative adversarial network with feature matching and task-specific perceptual loss[END_REF], the experts were asked to perform two tasks. First, they assigned a score [START_REF] Malm | Simulation of bright-field microscopy images depicting Pap-smear specimen[END_REF][START_REF] Calimeri | Biomedical data augmentation using generative adversarial neural networks[END_REF][START_REF] Apou | Synthesizing whole slide images[END_REF][START_REF] Glotsos | Image quilting and wavelet fusion for creation of synthetic microscopy nuclei images[END_REF][START_REF] Vitale | Improving realism in patientspecific abdominal ultrasound simulation using cyclegans[END_REF] for each synthesized image, where range 1-3 was considered to be low quality and 4-5 high quality. Second, the experts gave the amyloid status (positive vs. negative) for each image, the amyloid status being a biomarker used in the differential diagnosis of dementia. The status defined on the standard-dose ground truths and the synthesized images were compared. The consistency between the amyloid status showed whether the method was able to maintain the pathological features.

As stated in [START_REF] Xu | An empirical study on evaluation metrics of generative adversarial networks[END_REF], several works try to convince the readers about correlation of proposed metrics with human evaluation [START_REF] Salimans | Improved techniques for training gans[END_REF][START_REF] Lopez-Paz | Revisiting classifier two-sample tests[END_REF]. However, they also state that expert evaluation can be biased towards the visual quality of synthesized images and neglect the overall distributional characteristics, which are important for unsupervised learning.

Pairwise Comparison

In medical image synthesis, a majority of approaches require paired images in their training process, for example when learning to synthesise CT from MR images [START_REF] Burgos | Attenuation correction synthesis for hybrid PET-MR scanners: Application to brain studies[END_REF][START_REF] Dowling | Automatic Substitute Computed Tomography Generation and Contouring for Magnetic Resonance Imaging (MRI)-Alone External Beam Radiation Therapy From Standard MRI Sequences[END_REF][START_REF] Han | MR-based synthetic CT generation using a deep convolutional neural network method[END_REF][START_REF] Wolterink | Deep MR to CT synthesis using unpaired data[END_REF][START_REF] Nie | Medical Image Synthesis with Context-Aware Generative Adversarial Networks[END_REF], generate MR images of a certain sequence from MR images of another sequence [START_REF] Roy | Magnetic Resonance Image Example-Based Contrast Synthesis[END_REF][START_REF] Jog | Random forest regression for magnetic resonance image synthesis[END_REF][START_REF] Dar | Image Synthesis in Multi-Contrast MRI With Conditional Generative Adversarial Networks[END_REF], denoise low-dose CT images [START_REF] Wolterink | Generative Adversarial Networks for Noise Reduction in Low-Dose CT[END_REF][START_REF] Yang | Low-Dose CT Image Denoising Using a Generative Adversarial Network With Wasserstein Distance and Perceptual Loss[END_REF] or perform super-resolution [START_REF] Chen | Efficient and Accurate MRI Super-Resolution Using a Generative Adversarial Network and 3D Multi-level Densely Connected Network[END_REF]. Cross-modality synthesis is also present in microscopy imaging, where the attempt to reduce time-consuming and laborious tissue preparation results in synthesizing fluorescence images from the bright-field pairs [START_REF] Lee | Deephcs: Brightfield to fluorescence microscopy image conversion using deep learning for label-free high-content screening[END_REF]. A consequence of the need for paired images during training is that both reference and synthetic images are also often available for evaluation. The quality of the pairwise estimates is typically controlled by measuring the difference to the so called reference image (also known as ground truth or annotation) in the pair. The term reference image is also a reason why the pairwise comparison is sometimes called a full reference image quality assessment.

Generic pairwise performance measures

The measures the most frequently used to evaluate the synthesis accuracy by comparing real and synthetic images in a pairwise manner are listed below.

Mean absolute error

One of the most common measures is the mean absolute error (MAE) [START_REF] Burgos | Attenuation correction synthesis for hybrid PET-MR scanners: Application to brain studies[END_REF][START_REF] Dowling | Automatic Substitute Computed Tomography Generation and Contouring for Magnetic Resonance Imaging (MRI)-Alone External Beam Radiation Therapy From Standard MRI Sequences[END_REF][START_REF] Han | MR-based synthetic CT generation using a deep convolutional neural network method[END_REF][START_REF] Wolterink | Deep MR to CT synthesis using unpaired data[END_REF][START_REF] Nie | Medical Image Synthesis with Context-Aware Generative Adversarial Networks[END_REF][START_REF] Hiasa | Cross-modality image synthesis from unpaired data using CycleGAN -effects of gradient consistency loss and training data size[END_REF]. It is defined as the absolute difference between intensities in pixels of the simulated and ground truth image:

M AE = n i=1 |y i -x i | n , (3)
where x i and y i are the intensity values of the i-th pixel/voxel of the real and synthesized image, respectively, and n is the number of pixel/voxel pairs. In the context of PET attenuation correction [START_REF] Burgos | Attenuation correction synthesis for hybrid PET-MR scanners: Application to brain studies[END_REF] or MR-only radiotherapy treatment planning [START_REF] Dowling | Automatic Substitute Computed Tomography Generation and Contouring for Magnetic Resonance Imaging (MRI)-Alone External Beam Radiation Therapy From Standard MRI Sequences[END_REF], the MAE is often chosen as error metric as it is well suited when comparing CT images due to their quantitative nature.

Peak signal-to-noise ratio

The peak signal-to-noise ratio (PSNR) is a measure derived from the mean squared error (MSE):

M SE = n i=1 (y i -x i ) 2 n , (4)
with the same notations as for the MAE: x i and y i are the intensity values of the i-th pixel/voxel of the real and synthesized image, respectively, and n is the number of pixel/voxel pairs. PSNR compares the maximum of the intensity in the image with the error between the estimated and the ground truth image given by the MSE:

P SN R = 10 log 10 M AX 2 M SE , (5)
where M AX is the maximum possible intensity of the image. Higher values of PSNR relate to better simulation. The application of PSNR to validate the plausibility of synthetic images is apparent in many papers [START_REF] Wolterink | Deep MR to CT synthesis using unpaired data[END_REF][START_REF] Nie | Medical Image Synthesis with Context-Aware Generative Adversarial Networks[END_REF][START_REF] Jog | Random forest regression for magnetic resonance image synthesis[END_REF][START_REF] Dar | Image Synthesis in Multi-Contrast MRI With Conditional Generative Adversarial Networks[END_REF][START_REF] Wolterink | Generative Adversarial Networks for Noise Reduction in Low-Dose CT[END_REF][START_REF] Yang | Low-Dose CT Image Denoising Using a Generative Adversarial Network With Wasserstein Distance and Perceptual Loss[END_REF][START_REF] Chen | Efficient and Accurate MRI Super-Resolution Using a Generative Adversarial Network and 3D Multi-level Densely Connected Network[END_REF][START_REF] Lee | Deephcs: Brightfield to fluorescence microscopy image conversion using deep learning for label-free high-content screening[END_REF][START_REF] Hiasa | Cross-modality image synthesis from unpaired data using CycleGAN -effects of gradient consistency loss and training data size[END_REF].

Structural similarity

Even though they are simple to compute, the ability of the MAE, MSE and PSNR measures to perceive visual quality is limited [START_REF] Wang | Image quality assessment: From error visibility to structural similarity[END_REF]. Exploiting known characteristics of the human visual system, Wang et al. [START_REF] Wang | Image quality assessment: From error visibility to structural similarity[END_REF] proposed the structural similarity (SSIM), which compares local patterns of pixel intensities that have been normalised for luminance and contrast. SSIM has been widely adopted by the image synthesis community [START_REF] Jog | Random forest regression for magnetic resonance image synthesis[END_REF][START_REF] Dar | Image Synthesis in Multi-Contrast MRI With Conditional Generative Adversarial Networks[END_REF][START_REF] Yang | Low-Dose CT Image Denoising Using a Generative Adversarial Network With Wasserstein Distance and Perceptual Loss[END_REF][START_REF] Chen | Efficient and Accurate MRI Super-Resolution Using a Generative Adversarial Network and 3D Multi-level Densely Connected Network[END_REF][START_REF] Lee | Deephcs: Brightfield to fluorescence microscopy image conversion using deep learning for label-free high-content screening[END_REF]. SSIM is computed as a function of three components, luminance, contrast and structure, as follows:

SSIM (x, y) = (2µ x µ y + C 1 )(2σ xy + C 2 ) (µ 2 x + µ 2 y + C 1 )(σ 2 x + σ 2 y + C 2 ) , (6)
where x is the simulated image, y is the ground truth image, µ i is the mean value of image i, σ i is the variance of image i, σ xy is the covariance of images x and y. The constants C j are included to avoid instability when µ 2 x + µ 2 y , respectively σ 2

x + σ 2 y , is very close to zero. Therefore, they are set to C j = (K j • L) 2 , where L is the dynamic range of the pixel values (e.g. 255 for 8-bit grayscale images), and K ≪ 1 is a small constant. The original values of K j were set to K 1 = 0.01 and K 2 = 0.03 [START_REF] Wang | Image quality assessment: From error visibility to structural similarity[END_REF]. The higher the value of SSIM, the better the quality of the synthesized image.

In [START_REF] Hore | Image quality metrics: Psnr vs. ssim[END_REF], the relationship between SSIM and PSNR was investigated with the conclusion that the values of the PSNR can be predicted from the SSIM and vice-versa. Additionally, the PSNR and SSIM mainly differ in the sensitivity to image degradations. Similar conclusions can be found in [START_REF] Dosselmann | A comprehensive assessment of the structural similarity index[END_REF], where the association between SSIM and MSE was shown.

In their work, Mason et al. [START_REF] Mason | Comparison of Objective Image Quality Metrics to Expert Radiologists' Scoring of Diagnostic Quality of MR Images[END_REF] measured the correlations between ten pairwise metrics and the subjective score of five radiologists when assessing the quality of MR images. They showed that metrics such as SSIM and PSNR are potentially not ideal surrogate measures of MR image quality as determined by radiologist evaluation.

At this point, the reader can seriously wonder whether it is manageable to perform a proper data validation with all the previously-described metrics, together with all those that will be introduced in the following text. Fortunately, it is not needed. The above mentioned metrics (MAE, PSNR, SSIM) are the most common distance metrics and are accepted as a de-facto standard among the pair-wise metrics. Therefore, if you plan to generate some dataset that Biomedical Image Synthesis and Simulation, Chapter 25 form pairs (e.g. MRI-CT), you will likely be expected to evaluate at least one of these. All the metrics that come in the following paragraphs are not less important but rather mostly proposed/derived for specific purposes, so their usage is somewhat limited.

Application-specific pairwise performance measures

Pseudo-healthy image synthesis

To detect anomalies and better understand changes induced by diseases, Xia et al. [START_REF] Xia | Pseudo-healthy synthesis with pathology disentanglement and adversarial learning[END_REF] proposed to create subject-specific pseudo-healthy images from pathological ones using a CycleGAN. They assessed the quality of the image synthesis process using generic metrics, but they also designed new metrics tailored to their application: the so called healthiness, identity, and deformation correction metrics.

Healthiness The healthiness (h) expresses a fraction of (unwanted) pathology areas present in pseudo-healthy images and can be computed as:

h = 1 - E xp∼P [N (f p (G(x p )))] E mp∼Pm [N (f p (x p ))] (7)
where x p is the pathological image, f p the function providing a segmentation of pathological regions, G(•) the CycleGAN deriving the pseudo-healthy image from the pathological image, m p the ground truth mask of pathological regions and N (•) the number of pixels labelled as pathological by f p . One should pay attention that the healthiness can be strongly influenced by the quality of segmentation function f p .

Identity To measure the CycleGAN ability to preserve the subject identity (iD), one can measure the structural similarity of the original real pathology image and the derived pseudo-healthy image outside the pathological regions:

iD = MS-SSIM[(1 -m p ) ⊙ G(x p ), (1 -m p ) ⊙ x p ] (8)
where MS-SSIM stands for multiscale structural similarity [START_REF] Wang | Multi-scale structural similarity for image quality assessment[END_REF], x p is a pathological image, m p is its corresponding pathology mask, G(•) is the generator of pseudo-healthy images, and ⊙ is the pixel-wise multiplication.

Deformation correction In the case where brain tissue has recovered after some surgery or noninvasive therapy, there may be apparent structural changes. The deformation correction measure aims to assess whether such deformations have taken place. To avoid brightness influence, the images are first converted into edge maps. The classifier, that was trained over the set of Biomedical Image Synthesis and Simulation, Chapter 25 edge maps of healthy images, is used to judge whether the image of treated brain contains some deformations. The output of such classifier is a continuous number between 0 and 1 [START_REF] Xia | Pseudo-healthy synthesis with pathology disentanglement and adversarial learning[END_REF].

In this section, the validation methods used for pairs of images were reported. Note that the frameworks for synthesis of paired data are often based on algorithms using a loss function. The loss function should differ from the validation method since we want to avoid overfitting.

Dataset comparison

Some synthetic data are the result of frameworks generating a whole set of images according to the given input parameters. A set of generated images resulting from these tools should be comparable with the set of real images corresponding to the particular application of the synthesized data.

If we want to compare the images not only in the corresponding pixels but also in some quantitative characteristics, the use of image descriptors is suggested. Those image descriptors can represent characteristics such as color, shape, texture or some features from the frequency domain. A survey of image feature descriptors was published in [START_REF] Majtner | Texture-based image description in fluorescence microscopy[END_REF]. The values of a particular hidden layer of a neural network can also be considered as a descriptor (e.g. for the computation of the Fréchet Inception Distance -discussed further in 4.1.7). The most commonly derived descriptors used in validation methods of synthetic images are:

• Haralick texture descriptors

Haralick descriptors [START_REF] Haralick | Textural features for image classification[END_REF] represent a set of texture descriptors derived from so called gray level co-occurrence matrix. In [START_REF] Svoboda | MitoGen: A Framework for Generating 3D Synthetic Time-Lapse Sequences of Cell Populations in Fluorescence Microscopy[END_REF], contrast, correlation, homogeneity, and maximal correlation coefficient are reported to show the similarity between the real and synthetic image data. Haralick descriptors are an input parameter to affinity propagation in [START_REF] Kovacheva | A model of the spatial tumour heterogeneity in colorectal adenocarcinoma tissue[END_REF].

• Central moments A number of central moments, e.g. variance, skewness, kurtosis, can be calculated over pixel intensities. They were used for instance in [START_REF] Svoboda | On simulating 3D fluorescent microscope images[END_REF][START_REF] Svoboda | Generation of 3D digital phantoms of colon tissue[END_REF][START_REF] Malm | Simulation of bright-field microscopy images depicting Pap-smear specimen[END_REF].

• Subcellular location features [START_REF] Boland | A neural network classifier capable of recognizing the patterns of all major subcellular structures in fluorescence microscope images of HeLa cells[END_REF] This respectably extensive pack of descriptors include various types of patterns, such as texture features, Zernike moment features, object skeleton features and many others. In [START_REF] Zhao | Automated learning of generative models for subcellular location: Building blocks for systems biology[END_REF], the contribution of features on the classification of multiple cell objects was compared.

Biomedical Image Synthesis and Simulation, Chapter 25

• Local binary pattern [START_REF] Ojala | Performance evaluation of texture measures with classification based on Kullback discrimination of distributions[END_REF] Local binary pattern (LBP) is a texture operator which assigns to each pixel of an image a binary 8-bit number by thresholding its eight neighbors. LBP has been applied to ultrasound images to validate the method described in [START_REF] Singh | Texture based quality analysis of simulated synthetic ultrasound images using Local Binary Patterns[END_REF].

• Mean square displacement [START_REF] Michalet | Mean square displacement analysis of single-particle trajectories with localization error: Brownian motion in an isotropic medium[END_REF] The mean square displacement (MSD) descriptor is used for the validation of tracking objects. Svoboda et al. [START_REF] Svoboda | MitoGen: A Framework for Generating 3D Synthetic Time-Lapse Sequences of Cell Populations in Fluorescence Microscopy[END_REF] compared MSD using a histogram, displacement profiles and ensemble-average MSD curves.

• Spectral and spatial sharpness [START_REF] Vu | S3: A spectral and spatial sharpness measure[END_REF] In [START_REF] Zhao | Smore: A self-supervised anti-aliasing and super-resolution algorithm for mri using deep learning[END_REF], the spectral and spatial sharpness accompanied the standard SSIM and PSNR metrics to validate the quality of super-resolved MR images.

The above-mentioned list of descriptors is not limited. One can easily derive a data-specific descriptor that fits particular needs. For example, in the field of live cell analysis, the derived features such as the number of tracks, the average number of time points in tracks, or the average speed of particles may help to effectively compare the time-lapse image sequences [START_REF] Paavolainen | Application independent greedy particle tracking method for 3D fluorescence microscopy image series[END_REF].

One should keep in mind that image descriptors need not be necessarily used for the comparison of synthetic and real datasets. In case one is missing the real data, it is common to analyze solely the synthetic data and study the behaviour of the measured characteristics. This approach is also known as no reference image quality assessment.

The image descriptors are commonly used as input parameters to the simulation frameworks (number of elements/objects, shape, speed of motion, etc.) or are an integral part of loss functions in GANs or variational autoencoders. One should pay attention that the use of image descriptors as the inputs for the simulation together with the subsequent validation of the same features over the generated data is not only meaningful but may be even misleading. The results would be biased and therefore worthless.

Below, the reader can find the detailed explanation of the most prominent methods that utilize image descriptors and that are relevant to the comparison of sets of biomedical images.

Measures based on a comparison of probability distributions

Kolmogorov-Smirnov test

The Kolmogorov-Smirnov test (K-S test) [START_REF] Massey | The Kolmogorov-Smirnov test for goodness of fit[END_REF] is used when comparing the empirical distribution of a quantitative variable with a particular theoretical Biomedical Image Synthesis and Simulation, Chapter 25

distribution (e.g. the normal distribution). The K-S test can be used as well when comparing two samples and their empirical distributions (not necessarily lying under certain model distributions) [START_REF] Sorokin | FiloGen: A Model-Based Generator of Synthetic 3D Time-Lapse Sequences of Single Motile Cells with Growing and Biomedical Image Synthesis and Simulation, Chapter 25 Branching Filopodia[END_REF][START_REF] Svoboda | MitoGen: A Framework for Generating 3D Synthetic Time-Lapse Sequences of Cell Populations in Fluorescence Microscopy[END_REF].

The tested null hypothesis is that the cumulative distribution functions are the same for both samples A = {a i |i = 1, . . . , n A ; a i ∈ R} and B = {b i |i = 1, . . . , n B ; b i ∈ R}. The empirical cumulative distribution for sample A is defined as

F A (x) = #{a∈A|a≤x} n A
and in the same manner for sample B. In our case, A can stand for a sample of real data and B for a sample of synthetic data. The test statistic for the two samples is based on the largest distance (see Fig. 2) between the two empirical distribution functions, which is

KS(A, B) = sup x∈R |F A (x) -F B (x)| . ( 9
)
The value of the Kolmogorov-Smirnov statistic is then compared to the critical value and the null hypothesis is rejected if KS(A, B) >

1 n A + 1 n B κ α
, where κ α is chosen according to the level of significance α. Note that the Kolmogorov-Smirnov test is a non-parametric technique, which means it has no assumption put on the given data.

Kullback-Leibler divergence

The Kullback-Leibler (K-L) divergence, sometimes called relative entropy [START_REF] Do | Texture similarity measurement using Kullback-Leibler distance on wavelet subbands[END_REF], belongs to the family of probability measures. It is a measure of the dissim-Biomedical Image Synthesis and Simulation, Chapter 25 ilarity between two random quantities, in particular between two probability measures, as the Kolmogorov-Smirnov distance [START_REF] Venturini | Statistical distances and probability metrics for multivariate data, ensembles and probability distributions[END_REF].

The K-L divergence measures the inefficiency of assuming that a given distribution is P B (x) when the true distribution is P A (x). It is a kind of Bregman divergence and it is defined for discrete samples as:

D KL (A ∥ B) = x∈X P A (x) log P A (x) P B (x) , (10)
where P A (x) and P B (x) are two probability distributions [START_REF] Mackay | Information Theory, Inference and Learning Algorithms[END_REF].

Compared with statistical distances, the statistical divergence is not symmetric, i.e., it is not a metric. The K-L divergence returns the value of divergence of two probability distributions. It reveals 0 if and only if P A (x) = P B (x) otherwise it is non-negative.

This measure was applied and visualized using boxplots for comparison of synthetic and real data of microscopy nuclei images in [START_REF] Glotsos | Image quilting and wavelet fusion for creation of synthetic microscopy nuclei images[END_REF].

Kernel maximum mean discrepancy

Maximum mean discrepancy (MMD) [START_REF] Gretton | A kernel method for the two-sample-problem[END_REF] is another statistic measuring the difference between two probability distributions on the basis of samples drawn from each of them. More specifically, the aim is to find smooth functions resembling the sample values.

The test statistic is the difference between the mean function values in terms of kernel functions:

M M D 2 (P r , P q ) = E xr,x ′ r ∼Pr xq,x ′ q ∼Pq k(x r , x ′ r ) -2k(x r , x q ) + k(x q , x ′ q ) (11)
where P r and P q are the probability distributions for some fixed kernel function k, where r and r ′ are independent variables with distribution P r ; q and q ′ are independent variables with distribution P q . The biased empirical estimate of MMD is given by substitution of empirical estimates for expected values to:

M M D[R, S] = 1 m 2 m i,j=1 k(x i , x j ) - 2 mn m,n i,j=1 k(x i , y j ) + 1 n 2 n i,j=1 k(y i , y j ) 1 2
(12) where R and S (|R| = m, |S| = n) are two samples to be compared.

Large differences indicate that samples are from different distributions, whereas small MMD suggests that the distributions are identical. In their empirical study, Xu et al. [START_REF] Xu | An empirical study on evaluation metrics of generative adversarial networks[END_REF] showed that MMD satisfies most of the desirable properties and provided that the distances between samples are computed in a suitable feature space. To reach this conclusion, they designed a score that includes dropping and collapsing modes and is able to detect overfitting.

MMD was applied for quality assessment of generated data for example in [START_REF] Kang | Visual and quantitative evaluation of amyloid brain pet image synthesis with generative adversarial network[END_REF] and [START_REF] Kwon | Generation of 3d brain mri using autoencoding generative adversarial networks[END_REF].

Mutual information

Mutual information (MI) is related to joint entropy as introduced simultaneously by Viola [START_REF] Viola | Alignment by maximization of mutual information[END_REF] and Maes [START_REF] Maes | Multimodality image registration by maximization of mutual information[END_REF].

Given the images I and J, the joint entropy is defined as

H(I, J) = - a,b p IJ (a, b) log p IJ (a, b) (13)
where p IJ is the joint probability distribution of pixel intensities associated with images I and J. The joint entropy is minimized when the images I and J are similar. The individual (marginal) entropy for I is

H(I) = - a p I (a) log p I (a) (14)
and for J in the same manner. MI is then given by

M I = H(I) + H(J) -H(I, J) . (15)
In [START_REF] Hiasa | Cross-modality image synthesis from unpaired data using CycleGAN -effects of gradient consistency loss and training data size[END_REF], mutual information was utilized to measure the similarity between original MR image and synthesized CT image resulting from a CycleGAN [START_REF] Zhu | Unpaired image-to-image translation using cycle-consistent adversarial networks[END_REF], when the paired reference image was not available.

Regional mutual information

Regional mutual information (RMI) was introduced by [START_REF] Russakoff | Image similarity using Mutual Information of Regions[END_REF]. Compared with mutual information, the RMI takes into account information from the neighborhood of each pixel.

RMI is computed for two images I and J (each with corresponding pixels [I ij , J ij ]) in the following manner:

1. For each of the pixels [i, j], a vector v ij is created with the values coming from the two co-occurence matrices [START_REF] Haralick | Textural features for image classification[END_REF] from a neighborhood of radius r. The vector v ij is then considered as a point p k in a d-dimensional space where d = (2r + 1) 2 . Given the radius r and m × n images, there is a distribution of N = (m -2r)(n -2r) points represented by a d × N matrix P = [p 1 , ..., p N ].

2. The points are consequently centered by the mean:

P 0 = P -1 N N k p k . 3.
The covariance of the points is then computed as C = 1 N P 0 P T 0 .

Biomedical Image Synthesis and Simulation, Chapter 25

RM I = H g (C I ) + H g (C J ) + H g (C) . ( 16
)
As it is also mentioned in [START_REF] Russakoff | Image similarity using Mutual Information of Regions[END_REF], the corresponding pixels in the edges of the image can be handled in a number of ways.

In [START_REF] Glotsos | Image quilting and wavelet fusion for creation of synthetic microscopy nuclei images[END_REF], the authors used RMI pair-wise, not only for the real and synthetic counterparts, but also for each pair of real-synthetic images, and also for pairs of real images to express the measure of "alikeness". For the vizualization of the results, the boxplots were used for RMI values between real images next to boxplot of RMI between real and synthetic images.

Inception score

The inception score (IS) [START_REF] Salimans | Improved techniques for training gans[END_REF] is a measure for objective evaluation of trained GANs which was designed to correlate very well with subjective human judgement.

A deep convolutional network model (Google Inception v3 network [START_REF] Szegedy | Rethinking the inception architecture for computer vision[END_REF]), pre-trained on the large scale ImageNet dataset [START_REF] Deng | Imagenet: A large-scale hierarchical image database[END_REF], is used to classify the generated images from the proposed trained model G. The probabilities of images belonging to particular classes are computed and used for the evaluation of IS. The IS for a trained model G is computed as:

IS(G) = exp E X∼pg KL(p(y|X) ∥ p(y)) (17)
where X ∼ p g means that an image X is sampled from p g , p(y|X) is the conditional label distribution (class label y conditional on image X), and p(y) = X p(y|X)p g (X) is the marginal class distribution. KL is the K-L divergence (as defined in 4.1.2), E is an expected value, which is further exponentiated for an easier comparison. Two main properties are covered in IS -image quality (similarity of an image to a specific object) and image diversity (whether the dataset contains a wide range of generated objects). The distribution should have a low entropy for meaningful objects. Also, the marginal integral p(y|X = G(z))dz should have high entropy to provide varied images by the generator. The lowest value of IS is 1, the highest value is the number of classes defined in the pre-trained model. The calculation of IS assumes a large enough number of samples (i.e. 50000 images) to be able measure the diversity.

As a validation method, IS was applied in [START_REF] Wang | Semi-supervised mp-mri data synthesis with stitchlayer and auxiliary distance maximization[END_REF] for two proposed models designed to generate multi-parameter MRI data. IS was also evaluated in [START_REF] Calimeri | Biomedical data augmentation using generative adversarial neural networks[END_REF] to show the quality of synthesized MRI brain slices using GAN.

Fréchet distance and Fréchet inception distance

The Fréchet distance [START_REF] Fréchet | Sur quelques points du calcul fonctionnel[END_REF] is a measure of distance between curves that takes into account the location and ordering of the points along the curves. It is well known for its simile to walking a dog on a leash, who is crossing the curved path. A discrete variant of this measure, which can be evaluated over any feature vectors, exists [START_REF] Eiter | Computing discrete frechet distance[END_REF]. However, the most popular application in the field of image synthesis is the use of Inception feature vectors from Inception v3 model (called as the Fréchet inception distance (FID)), performed for example in [START_REF] Kang | Visual and quantitative evaluation of amyloid brain pet image synthesis with generative adversarial network[END_REF][START_REF] Xu | Bpgan: Bidirectional ct-to-mri prediction using multi-generative multi-adversarial nets with spectral normalization and localization[END_REF]. It was applied also on SonoNet in [START_REF] Lee | Generating controllable ultrasound images of the fetal head[END_REF] as Fréchet SonoNet distance.

FID summarizes the distance between the Inception feature vectors for real and generated images in the same domain. More specifically, the pretrained classification model is applied to both the real and generated images and the Inception feature vectors from the hidden layer are extracted, thus the final classification probabilities are not used. It was proposed in [START_REF] Heusel | Gans trained by a two time-scale update rule converge to a local nash equilibrium[END_REF] as an improvement over the existing IS (described in section 4.1.6), which also uses the Inception feature vectors v3 for evaluating the quality of generated images by GAN models. They showed in the paper that FID is more consistent with the noise level than IS. FID is a distance while IS is a score with a maximum value.

The definition of FID is given by:

d 2 = ||µ 1 -µ 2 || 2 + T r(C 1 + C 2 -2 * C 1 * C 2 ), ( 18
)
where µ 1 and µ 2 refer to the feature-wise mean of the real and generated images, where each element is the mean feature observed across the images. C 1 and C 2 are the covariance matrices for the real and generated feature vectors.

T r is a trace -a linear algebra operation -the sum of the elements on a diagonal of the square matrix. The best value of FID is 0.0 indicating that the two groups of images are identical, thus low FID values mean that the two groups of images are similar (or they have similar statistics).

Measures based on clustering and classification

Affinity propagation

Affinity propagation (AP) [START_REF] Frey | Clustering by passing messages between data points[END_REF] is a method of clustering data. The data points create a network of nodes, where each data point is a potential exemplar. The data points given to the input of this procedure are iteratively examined according to the measure of similarity (e.g. squared error in Euclidean distance) until a good set of exemplars and corresponding clusters emerges.

Initially, the input matrix s is computed over all combinations and each "similarity" is set to some optimization criterion, for example a negative squared Biomedical Image Synthesis and Simulation, Chapter 25 error:

s(i, k) = -∥x i -x k ∥ 2 (19)
for points x i and x k . The values of elements belonging to the diagonal (s(i, i)) are able to prefer some of the points to be an exemplar. Two types of "messages" are evaluated in each step for two points: the "responsibility" message r(i, k) sent from point i to point k and the "availability" a(i, k) sent from point k to point i. The first message, "responsibility":

r(i, k) ← s(i, k) -max k ′ ̸ =k {a(i, k ′ ) + s(i, k ′ )} (20)
reflects the actual evidence for how well-suited point k is to serve as the exemplar for point i (regarding the potential exemplars of point i). The second message, "availability"

a(i, k) ← min    0, r(k, k) + i ′ / ∈{i,k} max{0, r(i ′ , k)}    (21)
reflects the accumulated evidence for how appropriate it would be for point i to choose point k as its exemplar (regarding the support from other points that point k should be an exemplar). These are initially set to zeros. The "self-availability" a(k, k) is updated in each step as

a(k, k) ← i ′ ̸ =k max{0, r(i ′ , k)} . (22)
The messages are sent until either the cluster boundaries remain unchanged over a number of iterations, or after some predetermined number of iterations. The advantages of AP are that no specific number of clusters has to be predefined.

In [START_REF] Kovacheva | A model of the spatial tumour heterogeneity in colorectal adenocarcinoma tissue[END_REF], affinity propagation was applied to [START_REF] Lehmussola | Computational framework for simulating fluorescence microscope images with cell populations[END_REF] Haralick texture features of 20 images of real nuclei to get different phenotypes of images. These were compared to the phenotypes clustered on synthetic images. The frequencies of phenotypes in both groups were compared using histograms.

Classification

We should keep in mind that the aim of validation methods is not to find the differences between the synthetic and real data, but to explore the similarity of the synthetic data compared with the real data. One could suggest that a classification task distinguishing real and synthetic data is a suitable solution. However, the bad result of any classifier does not prove that the data are not similar. For example in Fig. 3, it is not possible to distinguish the two classes of points (blue and yellow) with any discrimination function, but at first sight one can see that the groups are not homogeneous. Furthermore, the result of a classification could be confused with a bad choice of classifier. Some works employed classification not for differentiation of real versus synthetic dataset, but for verification that the synthetic images will be classified to the appropriate subcategories as well as real data. For instance, in [START_REF] Zhao | Automated learning of generative models for subcellular location: Building blocks for systems biology[END_REF] the synthetic images for ten classes of organelles were generated. The support vector machine classification was applied to both real and synthetic data, and results (such as confusion matrix) were reported for comparison. Similarly, in [START_REF] Lehmussola | Synthetic images of high-throughput microscopy for validation of image analysis methods[END_REF] the k-nearest neighbor to classify bacteria classes based on their shape was performed.

Measures based on a visualization of transformed data

Histogram

A histogram is a widely used plot for visualization of absolute or relative frequencies of a categorized quantitative variable. The quantitative variable (e.g., an image descriptor) determines the axis x with a specific binning (categorization), and the axis y counts the frequency -relative (f ) or absolute (f * ) -in the bins:

f (j) = n j /n d j , f * (j) = n j d j , ( 23
)
where n is the total number of observations, n j is the number of observations for a particular category j, and d j is the width of this category. The shape of Biomedical Image Synthesis and Simulation, Chapter 25

the histogram should approximate the probability density function. However, the approximation is very sensitive to the choice of the number of bins.

A quality assessment of the synthesized images at the image descriptor level was used in [START_REF] Svoboda | MitoGen: A Framework for Generating 3D Synthetic Time-Lapse Sequences of Cell Populations in Fluorescence Microscopy[END_REF], where the distributions of five Haralick texture descriptors [START_REF] Haralick | Textural features for image classification[END_REF] were compared using histograms and the Kolmogorov-Smirnov test. In [START_REF] Kovacheva | A model of the spatial tumour heterogeneity in colorectal adenocarcinoma tissue[END_REF], histograms were used to compare the length of minor axis and ratio between minor and major axes of healthy crypts in colon tissue in real and synthetic images. The Gamma distribution of these two measures was fitted to the real data and subsequently plotted into both histograms.

Quantile-quantile plot

The quantile-quantile (Q-Q) plot [START_REF] Wilk | Probability plotting methods for the analysis for the analysis of data[END_REF] is designed for the visualization of distribution comparisons. The quantiles of two quantitative variables are plotted in the x-y figure against each other (see Fig. 4). The quantile matches (and therefore identical distribution) reveals if the points lie along the straight line with the slope of 1. Distant points from this line indicate deviation from the same distribution. The method is non-parametric with no assumption put on the input data. In case of image descriptors, the Q-Q plot can help to compare the distribution of both groups of data in each descriptor separately, i.e. univariately [START_REF] Sorokin | FiloGen: A Model-Based Generator of Synthetic 3D Time-Lapse Sequences of Single Motile Cells with Growing and Biomedical Image Synthesis and Simulation, Chapter 25 Branching Filopodia[END_REF][START_REF] Malm | Simulation of bright-field microscopy images depicting Pap-smear specimen[END_REF][START_REF] Svoboda | Generation of digital phantoms of cell nuclei and simulation of image formation in 3D image cytometry[END_REF][START_REF] Zhao | Automated learning of generative models for subcellular location: Building blocks for systems biology[END_REF].

Fourier ring correlation

The Fourier ring correlation (FRC) [START_REF] Nieuwenhuizen | Measuring image resolution in optical nanoscopy[END_REF] is a spatial frequency correlation function that measures the degree of correlation of two images at different spatial frequencies. The images are initially transformed to the frequency domain. Afterwards, the normalized average correlation is computed for N r concentric rings of increasing radius, which corresponds to increasing spatial frequencies centered around the (0,0) spatial frequency.

FRC for images I and J (which are transformed to the frequency domain) is computed as:

F RC(R) = i∈R I(r i ) • J(r i ) * ( i∈R |I(r i )| 2 )( i∈R |J(r i )| 2 ) ( 24
)
where R is the ring number and r i are the spatial coordinates. The values of the FRC draw a curve in a plot, which can be further investigated (see Fig. 5). This approach is found to be useful in electron and fluorescence microscopy [START_REF] Venkataramani | SuReSim: simulating localization microscopy experiments from ground truth models[END_REF][START_REF] Banterle | Fourier ring correlation as a resolution criterion for super-resolution microscopy[END_REF] for determining the resolution at which both images are consistent. Here, the 2σ criterion is specified as

F 2σ (R) = 2 N p (R)/2 (25)
where N p (R) is the number of pixels in the ring R. The optimal resolution is then found as the crossing of F RC and F 2σ curves. In [START_REF] Nečasová | Visual and quantitative comparison of real and simulated biomedical image data[END_REF], the validation method is based on the comparison of explained variability of both real and synthetic data in the same feature space. The descriptors (Haralick descriptors in this case) are initially preprocessed by principal component analysis to reduce the original number of dimensions into only three, easy-to-visualize, dimensions. The real and synthetic images are represented as data points in this feature subspace. Finally, the overlap of the clusters created around real and synthetic data (Fig. 6) is evaluated via Jaccard index as a quantitative measure of this technique. However this validation method assumes that the three principal components are able to explain the majority of the original feature space given by the descriptors.

t-distributed stochastic neighbor embedding

t-distributed stochastic neighbor embedding (t-SNE) [START_REF] Van Der Maaten | Visualizing high-dimensional data using t-sne[END_REF] is also a reduction (nonlinear in this case) from a multidimensional feature space into two or three dimensions designed for visualisation. The method is an improved variation of stochastic neighbor embedding (SNE) [START_REF] Hinton | Stochastic neighbor embedding[END_REF] with easier optimization and better results in terms of spread of the points in the 2D or 3D map. The images are represented as high-dimensional data points in a feature space. Like SNE, the Euclidean distances between data points are converted into conditional probabilities representing similarities between the images. The similarity of data point x j to data point x i in SNE is expressed as:

p j|i = exp(-∥x i -x j ∥ 2 /2σ 2 i ) k̸ =i exp(-∥x i -x k ∥ 2 /2σ 2 i ) (26)
where σ i is the variance of the Gaussian that is centered on data point x i . The probability corresponds to the probability that x j would be in a neighborhood of x i in proportion to their probability density under the Gaussian centered at x i . The values of p i|i are set to zero. The low-dimensional counterparts for data points x i , x j in a high-dimensional space are labeled as y i , y j . The conditional probability q j|i for the low-dimensional counterparts is evaluated in a similar way as in equation 26, but the variance of the Gaussian is set to 1 √ 2 . With the use of symmetric SNE, the joint probabilities are employed instead of conditional probabilities for a faster computation. They are set to p ij = p j|i +p i|j 2n

. Compared with SNE, t-SNE employs a heavy-tailed distribution -a Student t-distribution with one degree of freedom -instead of Gaussian in the low-dimensional map. The final joint probabilities q ij are defined as

q ij = (1 + ∥y i -y j ∥ 2 ) -1 k̸ =l (1 + ∥y k -y l ∥ 2 ) -1 . (27)
The particular locations of the data points in the map are determined by minimizing the K-L divergence over all data points using gradient descent. A validation of synthetic data using t-SNE was applied as a qualitative method in amyloid brain PET image synthesis by [START_REF] Kang | Visual and quantitative evaluation of amyloid brain pet image synthesis with generative adversarial network[END_REF]. Qualitative visual assessment of the real and synthetic distribution in the 2D plot via t-SNE was also applied in [START_REF] Han | Synthesizing diverse lung nodules wherever massively: 3d multi-conditional gan-based ct image augmentation for object detection[END_REF]. They compared two settings of training GAN and real images in the map. t-SNE was used in [START_REF] Sun | An adversarial learning approach to medical image synthesis for lesion detection[END_REF] to show the homogeneity of generated and real images with (or without) lesions. In this work, the aim was to generate normal-looking counterpart for the abnormal images with lesions. Similarly, in [START_REF] Xue | Selective synthetic augmentation with histogan for improved histopathology image classification[END_REF] t-SNE of the original and augmented cervical histopathology was a part of the qualitative validation of the training set. In [START_REF] Diaz-Pinto | Retinal image synthesis and semi-supervised learning for glaucoma assessment[END_REF], the authors used t-SNE in a different way -for a comparison among generated methods, not real images.

Measures applied to time-lapse sequences

Linear mixed models

Linear mixed models (LMM) also known as hierarchical linear models or linear models with mixed effects are an extension of linear models for dependent measurements. The extension is needed to overcome two important assumptions put on basic linear models, which are violated in data with dependent Biomedical Image Synthesis and Simulation, Chapter 25 measurements. One of them is that the residuals have to be of homogeneous variance and the other that the residuals should not be correlated. The dependency can result from longitudinal data, repeated measurements or clustered data, because one observation is related to another. This method is one of the possibilities for statistical comparison of time-lapse sequenced image data as well as other dependent datasets of images.

The linear mixed model [START_REF] West | Linear Mixed Models[END_REF] is considered to be of formula:

Y i = X i β + Z i b i + ε i , (28)
where Y i = (Y i1 , ..., Y in i ) ′ is the vector of repeated quantitative outcome measurements for object i, β is treated as a vector of fixed effects, i.e. populationaverage regression coefficients, and X i is a known design matrix linking β to Y i . The effect of β is the same for all the objects. In the case of validation of synthetic data, we are interested in the statistical significance of β j , which is one element of vector β. A coefficient β j is the estimated effect of X j , a binary variable holding the information about being a synthetic or a real image. Naturally, β j can be the only element of β if X j is the only explanatory variable in the design matrix X i . Random effects are held by b i , a vector of q object-specific regression coefficients. The columns in the Z i matrix represent observed values for the q predictor variables for the i-th subject, which have effects on the continuous response variable that vary randomly across subjects. In many cases, predictors with effects that vary randomly across subjects are represented in both the X i matrix and the Z i matrix.

The residuals ε i are distributed as ε i ∼ N (0, R i ), where R i is a covariance matrix and depends on i only through its dimension n i . The b i are distributed as b i ∼ N (0, D), independently of each other and of the ε i , where D is a covariance matrix of the random effects [START_REF] West | Linear Mixed Models[END_REF][START_REF] Laird | Random-effects models for longitudinal data[END_REF][START_REF] Molenberghs | A review on linear mixed models for longitudinal data, possibly subject to dropout[END_REF].

LMM were applied in [START_REF] Svoboda | Tubular network formation process using 3D cellular potts model[END_REF] to images of tubular network of endothelial cells. For assessment of the effect of being an image from a group of real or synthetic dataset, two standard measures were chosen: box counting fractal dimension and lacunarity [START_REF] Gould | Multifractal and lacunarity analysis of microvascular morphology and remodeling[END_REF][START_REF] Smith | Fractal methods and results in cellular morphology -dimensions, lacunarity and multifractals[END_REF], both describing the complexity of the structures depicted in the analyzed images.

Dynamic time-warping

Dynamic time warping (DTW) is an algorithm for measuring similarity between two time series which may vary in timing. In another words, one can differ from another only in being slower or faster in any part of the series ("warped"). The aim of the algorithm is to find an optimal alignment. The algorithm was originally applied for comparing speech patterns in automatic speech recognition [START_REF] Rabiner | Fundamentals of Speech Recognition[END_REF].

In [START_REF] Müller | Information Retrieval for Music and Motion[END_REF], the DTW for two time series X := (x 1 , x 2 , ...x N ) of length N ∈ N and Y := (y 1 , y 2 , ...y M ) of length M ∈ N is defined using a local cost measure.

The local cost measure is a function

c : F × F → R ≥0 , (29)
Finally, the optimal warping path is the one having minimal total cost among all possible warping paths. This task is a standard optimization problem which can be solved, e.g., by using dynamic programming. DTW as a validation method was applied in [START_REF] Svoboda | Image-based simulations of tubular network formation[END_REF]. On a dataset of tubular networks of epithelial cells, the descriptor computing the number of lagoons in the network was evaluated for each frame of the sequence. The curves resulted from real and synthetic image sequences showing the development of the networks were compared against each other.

Conclusion

A wide range of methods for the validation of generated data were described in this chapter. Even though the mathematical models offer a powerful tool in the inspection of synthetic data, they should be accompanied by a humandriven qualitative assessment. However, the human eye can be deceived by many factors, e.g. the level of details and number of dimensions, thus it is not sufficient as such. The qualitative and quantitative assessment should go hand-in-hand together and the method should be carefully chosen with regard to suitability for a particular application.

Figure 2 :

 2 Figure 2:The distance between two distribution functions for an equation of Kolmogorov-Smirnov test. In the one-sample case (left) the distance is measured between the specific theoretical cumulative distribution (black) and the empirical distribution (red). For the two-sample case (right) the distance is computed between the two empirical distribution functions.

Figure 3 :

 3 Figure 3: A classification problem. It is not possible to distinguish the two groups of points (blue and yellow) with any discriminator without transformation. However, at first sight, one can see that the groups are not homogeneous.

Figure 4 :

 4 Figure 4: Example of Q-Q plot in case of one-sample comparison (left) and twosample comparison (right). The optimal match of distributions is achieved when the points lie on the dashed line.

Figure 5 :

 5 Figure 5: Fourier ring correlation. The values of FRC for each ring are plotted as a curve (blue). The red curve depicts the F 2σ . Adopted from [76].

Figure 6 :

 6 Figure 6: Interactive visualizer showing real images (yellow) and synthetic images (blue) in the reduced feature space. The axes are the three main components given by principal component analysis. Each point corresponds to a particular image [77].

  comparing two features x, y from a fixed feature space F . Local cost measures for each pair of elements X and Y are composed into cost matrixC ∈ R N ×M defined by C(n, m) := c(x n , y m ) . (30)The goal is to find an alignment between X and Y with the minimal overall cost. A warping path is a sequence p = (p 1 , ...p L ) withp l = (n l , m l ) ∈ [1 : N ] × [1 : M ] for l ∈ [1 : L]that corresponds to the possibility to align X and Y with a particular total cost c p (X, Y ) with respect to the local cost measure c: c p (X, Y ) := L l=1 c(x n l , y m l ) .

Table 1 :

 1 Results of the experts' evaluation in[START_REF] Malm | Simulation of bright-field microscopy images depicting Pap-smear specimen[END_REF]. The classifiers considered included six experts and a random classifier. T P is the number of correctly classified real

images -true positives, T N is true negatives -correctly classified synthetic images, F P is false positives and F N is false negative cases.

  4. The joint entropy H g (C) is estimated. 5. The marginal entropies H g (C I ) and H g (C J ) are estimated where C I is a matrix in the top left of C and C J in the bottom right, both of size

	d 2 × d 2 .
	6. Finally, the RMI is calculated as

Biomedical Image Synthesis and Simulation, Chapter 25

Acknowledgments

We thank Yang Song from University of New South Wales for her invaluable feedback to this chapter.

N. Burgos received funding from the French government under management of Agence Nationale de la Recherche as part of the "Investissements d'avenir" programme reference ANR-10-IAIHU-0006 (Agence Nationale de la Recherche-10-IA Institut Hospitalo-Universitaire-6).