Coherence loss by sample dynamics and heterogeneity in x-ray laser diffraction

In its most common implementations, single particle imaging with x-ray free-electron lasers (XFELs) assumes the perfect coherence of the diffracted x-rays. However, XFEL imaging is destructive and exhibits a reduction of coherence due to sample dynamics (radiation damage). Sample heterogeneity also manifests as effective loss of coherence in averaged data, which is increasingly encountered in the development of 3D XFEL imaging techniques. Here we review coherence loss in XFEL imaging experiments, which has focused primarily on femtosecond electronic dynamics. We present extensions of coherence theory to XFEL crystallography, ion diffusion, the Coulomb explosion of a single particle and to heterogeneity using a spherical model. In each case, we provide illustrative simulations of the effects of coherence loss and discuss when these effects will need to be taken into account.


Introduction
Single particle imaging (SPI) with x-ray free-electron lasers (XFELs) is being actively developed for its potential applications to single virus particles and macromolecules [1]. These methods are based on coherent diffractive imaging (CDI) techniques [2] that measure far-field diffraction intensity, which is inverted to recover an image by retrieving the phase. The inversion is performed with algorithms that were originally developed for optical and synchrotron CDI [3,4]. In their simplest implementation, CDI algorithms use an assumption of full coherence. It has been recognized though that the destructive nature of the XFEL pulses induces electronic changes in the sample during the measurement, which reduces the coherence of the scattered x-rays [5,6].
Partial coherence has long been recognized as an important issue for synchrotron CDI [4]. The finite size of the effective source, the energy bandwidth and sample motion all reduce the coherence in synchrotron CDI experiments. CDI algorithms have been extended to account for the coherence properties of the beam [7] and multiple wavelengths [8,9]. It was shown that accounting for these effects significantly improves the robustness of the image reconstruction process. In the case of broadband CDI [9], the loss of coherence is compensated by the advantages of higher flux and lower exposure time. Sample dynamics on timescales faster than the duration of the measurement reduce the coherence of the scattered x-ray wave, even if the incident wave is fully coherent. Methods for imaging dynamic samples have been developed for multi-measurement scanning CDI methods, like ptychography [10,11], by exploiting the redundancy of the dataset. It is important to recognize that this does not provide a movie of a dynamical process, but rather images of the dominant structural correlations traversed by the sample during its motion.
In XFEL experiments sample dynamics are initiated by core-shell photoabsorption, which is approximately ten times more likely than elastic scattering [12] and starts from the first femtosecond of the pulse-sample interaction. This is followed by Auger decay and secondary ionization that occur within the first 10 fs [13]. Since ionization modifies individual atoms it has the greatest effect approaching atomic resolution. Ion motion by either Coulomb repulsion [14] or diffusion [15] only occurs after 10 fs and leads to a loss of shape and structure on nanometer and sub-nanometer length scales. Ion motion can potentially be outrun by sub-10-fs pulses. Damage depends on pulse parameters via the photoionization rate which is proportional to pulse intensity and decreases with photon energy [12]. The secondary processes, particularly secondary ionization, have a non-linear relationship to the photoionization rate [14] and numerical simulation is needed to model these effects.
Single XFEL pulses have been shown to be highly spatially and temporarily coherent [16,17]. As far as the incident pulse is concerned, the assumption of full coherence for CDI experiments is usually justified. However, sample-induced dynamics can reduce the coherence of the diffracted x-rays and it has been predicted that this will have a significant effect approaching atomic resolution [5]. There have been two studies of coherence loss due to electron dynamics induced by XFEL pulses [5,6]. These studies found that at current XFEL intensities a single coherent wave can be used to describe the scattered intensity, but that coherence effects are required to interpret the wave. In particular, the observed scattering strengths of different elements are modified by coherence effects. It was also shown that if the dynamical effects can be appropriately characterized [18], then they can be accounted for during the CDI analysis to recover an image of the undamaged structure [5].
Here we review the work that has been done so far on the coherence analysis of electron dynamics. We then proceed to consider new applications for coherence analysis. We show that coherence effects play a role in the apparent scattering strength of different elements in XFEL crystallography at high pulse intensities. We explore the use of coherence analysis to image the Coulomb explosion of a single molecule using a 1D spherical hydrodynamic model of the explosion.
It is also recognized that the biological samples of most interest are not perfectly reproducible and exhibit thermal variability or structural heterogeneity. Averaging data from multiple measurements of such samples is an unavoidable step in 3D imaging, but induce a loss of effective coherence due to the sample variability. The effects of thermal variability of a protein molecule on CDI has been studied with MD simulations [19], but the analysis did not include calculating effective coherence. Here, we consider the coherence analysis of data averaged from a heterogeneous sample and show that random variations of the size of spherical targets imposes an effective resolution limit on the image.

Background: coherence theory in diffractive imaging
The theory of coherence in diffractive imaging has been reviewed elsewhere [4]. Here we provide an introduction to the concepts of coherence theory that will be used in the rest of the article.
In optical coherence theory, a fully coherent x-ray beam is modeled by a complex scalar wave y r ( ). This is a monochromatic wave (temporally coherent) originating from a source that is effectively point-like (spatially coherent). The coherence of an x-ray beam is effectively reduced when it is composed of several waves, y a r; ( ), that vary with a parameter a and are not mutually coherent. The parameter a could, for example, denote wavelength, a point on the effective source or time (for a dynamical sample), and it is possible to generalize the formalism to multiple parameters. The beam can then be described by a mutual intensity function where w a ( ) is a weighting function that is required if y a r; ( ) is normalized. The theory of x-ray propagation and interaction with matter can be formulated using  ¢ r r , ; ( ) it is equal to the intensity of the beam when ¢ = r r. In practice, it is inconvenient to work with a six-dimensional form of  ¢ r r , ) can be diagonalised to generate what is known a modal representation: where the modes y a r ( ) are a set of orthogonal eigenfunctions and h a are real-valued weights. The size of the weights effectively determines the number of modes needed to described the beam, as weights that are close to zero correspond to modes that are not present in the beam. The modal representation is very effective for partially coherent beams with only a few modes, and it is the basis for extending coherent diffractive imaging algorithms to partially coherence beams.
In the context of XFEL imaging, the assumption of a fully coherent incident beam is justified for most practical purposes. There are two important situations where coherence can be reduced due to the sample. First there is the damage induced in the sample by the XFEL beam. In this case, the intensity is an incoherent average of the dynamically evolving sample: where q is the reciprocal space coordinate that is convenient for the description of far-field diffraction.
The second way the sample can reduce coherence occurs when measurements of a heterogeneous sample are averaged to improve signal to noise. This is a necessary step in single particle imaging, where large numbers of noisy two-dimensional snapshots are merged into a three-dimensional dataset. Since the measurement is destructive, each snapshot is taken from a fresh copy of the particle. If each of these particles is not identical (as assumed in the ideal case), then the merged three-dimensional intensity is an average of diffracted intensities over a set of sample structures  :

Coherence theory for electron dynamics and ion diffusion
There are two studies of coherence theory for XFEL-induced damage that have considered the effects of ionization [5,6].
Here we review this theory for single particle imaging and discuss the extensions to nuclear degrees of freedom. Coherence theory for single particle imaging begins with an atomic description of the scattering factor for the particle: Both the atomic scattering factors and the displacements are written as a function of time t to account for the effect of the XFEL pulse. The scattering factors are changed by the ionization and the resulting ion is displaced from its initial position by the Coulomb forces of surrounding ions. The parameter τ labels a damage scenario, which is a stochastic sequence of ionizations and ion displacements that will be different for each measurement. We can use t F t q, , ( ) to construct the mutual intensity matrix for XFEL damage, The average over damage scenarios, denoted by áñ t , is performed in practice by averaging many measurements, which is a necessary part of 3D imaging. Substituting equation (5) into equation (6) and rearranging leads to When written as a matrix, where N a is the number of atoms and N q is the number of q-space samples. The number of atoms in the smallest SPI targets are protein molecules with N a in the range of 10 10 3 5 -. The number of samples for a 3D Fourier space is, at minimum,´´= = N 100 100 100 10 . For a completely general damage scenario the matrix representation of is too large to diagonalize. Damage processes like ionization and ion diffusion which are not correlated between ions depend on the element of each ion, but not its initial position. For these processes, we can capture the essential physics by a matrix that depends only on the number of atomic species and the sampling of the magnitude of the scattering q. We assume the statistical independence of ionization and ion diffusion to write We further assume that for different ions, ¹ i j, ion diffusion is statistically independent. The ensemble average for the ion diffusion term can be evaluated to obtain In this case we can decompose the matrix repre- whereÃ is a smaller coherence matrix given by and T is a matrix with N N a q columns and N N z q rows that assigns each ion in the sample to an element and an initial position. Its elements are given by We can then diagonalizeÃ to obtain modes ) can be constructed from the modal representation ofÃ by writing This was the approach that was taken in earlier work [5] that modeled the time-dependence of t f t q, , Z ( ) due to photoabsorption. Interestingly, it was found that the diffracted x-rays could be characterized by a single dominant mode. The weighting of the second mode was 2-3% of the first mode. In principle this means that a CDI algorithm assuming full coherence can recover an image of the sample. However, electron dynamics affect the interpretation of this image and, in particular, the observed scattering strength of different elements. It was also shown that the dynamical effects can be corrected if the matrixÃ on a system of similar size and chemical composition, such a known protein structure.

Statistical simplifications for the coherence theory of induced electron-dynamics
The study of coherence and electron dynamics by Lorenz et al [6] used further simplifications based the statistical independence of ionization events on different atoms. When this independence holds, we have , , Z Z ( ) ( ) to simplify notation, we can use this result to define a matrix Ā by In the absence of ion diffusion, this matrix is related to (12), by is a Kronecker delta function. The advantage of this formulation is that the first term on the right hand side of equation (18) can be identified as q-dependent background scattering. Ā captures coherence properties of the remaining diffraction and is likely to be closer to a single mode representation thanÃ. Hence, this formulation suggests alternative ways to use information about electronic damage to improve the formulation of CDI algorithms by explicitly modeling a background component. Formulations based onÃ and Ā should be equivalent up to the validity of equation (16). Since the approximation in equation (16) is related to the approximation that underlies both approaches that was used to go from equation (8) to equation (12), we expect this not to separate the two models in practice.
The results of Lorenz et al broadly agree with those of Quiney et al in predicting one dominant mode for currently available experimental parameters. Notably, Lorenz et al used a rate-equations model to include secondary electron impact ionization from both escaping photoelectrons and electrons trapped in a thermalized bound-electron gas. It is a positive result for the future realization of CDI that these additional ionization mechanisms do not prevent a single mode formulation of single particle diffraction.

Coherence in XFEL crystallography
Radiation damage in XFEL crystallography is commonly studied by modeling the time-averaged atomic structure factors [15]. The coherent addition of many unit cells simultaneously averages over many damage scenarios in a single measurement, which may have led to the view that coherence analysis is not required. A more careful consideration shows this is not the case. Coherence loss in XFEL crystal diffraction is closely related to the theory outlined in section 4.
The scattering factor for the crystalline case can be written as where n is an index that ranges over all of the units cells in the crystal and N is the number of unit cells. The vector q hkl gives the location of a Bragg peak, which satisfies the condition for all values of n. Due to the stochastic nature of XFEL pulses by selfamplified spontaneous emission (SASE), the coherence time of XFEL pulses is less than a femtosecond [20]. As described in the introduction, most of the induced dynamical processes take of the order of a femtosecond or longer, which motivates constructing the mutual intensity function with an incoherent integration over pulse duration, as follows , , e , , e d .
In matrix form we can write this as where Ā is derived using the methods of section 4. We see that there is a coherence loss at the Bragg peaks due to the time integration over the duration of the pulse.
The coherence loss in crystals is very similar to coherence loss in single particles [5,6], so that the validity of the single mode approximation for single molecule imaging should carry over to the crystal case. In the simulations that follow, we have also found this to be true.
Even in a single mode approximation, the contributions of atomic species to the dominant mode are not necessarily equal to the average scattering factors [5]. Differences are expected to arise when atoms ionize at different rates. This point may not have been significant in previous studies because [5] focused primarily on photoionization, and in [6], valence excitation was considered with a rate-equation model, but only included lighter elements of similar atomic number (carbon, nitrogen and oxygen) that show similar ionization behavior. There are, however, crystallographic methods that use heavy elements explicitly for structure determination. For these heavier elements, coherence analysis will have more to say. We illustrate this by including sulfur in addition to the lighter elements H, C, N and O in a simulation that includes valence shell ionization. For sulfur we treat all sixteen electrons. The damage model we use is based originally on work by Hau-Riege [14] and is detailed in [21]. It includes photoionization, Auger decay, secondary ionization and ion diffusion. We simulate a 100 nm beam diameter and assume a sample volume of´100 100 100 nm 3 , which is assumed to contain the equivalent of 216 GroEL molecules. Figure 1 shows the differences for the scattering factor for carbon as function of pulse intensity, which indicates that as we approach and then exceed the maximum intensity currently achievable differences emerge between the average structure factor and the dominant coherent mode. Although only carbon is shown in figure 1, a similar trend is observed for other elements. If these differences were a common multiplicative factor, independent of atomic species and resolution, then it would not be significant in practice. This is not the case. For higher pulse intensities, there is variation as a function of q of up to 20%. At the highest intensity shown in  The first mode is weighted by the square root of its eigenvalue. The photon energy is 10 keV and the spot size is 100×100 nm 2 . The crystal size was´100 100 100 nm 3 and it was assumed to contain the equivalent of 216 GroEL molecules. figure 2, sulfurʼs contribution to the dominant mode scales more weakly than lighter elements. This result indicates that a more thorough treatment of coherence analysis in XFEL crystallography is needed for heavy elements such as sulfur, manganese, iron and other elements used in anomalous dispersion or isomorphic replacement techniques.

Coherence investigation of 1D radial expansion
The Coulomb explosion of a molecule in a high intensity XFEL beam sets in at around 10 fs. The explosion has been regarded as fatal for single particle imaging and it has been proposed to avoid it by using pulses less than 10 fs in duration. This limitation is unfortunate because current facilities produce pulses with optimal intensity in the range 30-50 fs. This limit may also be harsh, because it is has been shown that scattering of single particles from longer pulses contains a significant amount of information about the initial atomic positions [21]. The scattering intensity from an exploding sample is not consistent with the assumption of full coherence. Therefore, if it is possible to extend CDI methods to account for the effects of sample explosion, it must be done within the framework of partial coherence. Here we study the extent of the coherence loss in a 1D spherical model of the Coulomb explosion. This is a first step toward understanding if CDI can be extended to image an exploding sample.
The simulation of the 1D expansion divides the spherical model of the sample into radial layers. For each time-step it calculates the radial distribution of trapped electrons accounting for the Coulomb forces of both ions and electrons. For this calculation input is required from the rate-equation model, namely the total number trapped electrons and the distribution of ionization states for each element for each time-step. It then calculates the acceleration of the radial layers from the distributions of ions and trapped electrons and updates the radial position of the layers. The details of the simulation of the 1D radial expansion are given in [14,21]. As some of the photoelectrons escape the molecule, the total charge of the trapped electrons is less than the charge of the positive ions. The trapped-electron distribution tends to neutralize the core of the particle, leaving a outer shell of positively charged ions exposed to the Coulomb forces. In this way, the explosion proceeds by sequentially removing layers as shown in figure 3(a).
Once we have calculated the radial position of all the layers as a function of time, we can perform a modal analysis of the x-ray diffraction. The scattering factor for a sphere is the Fourier transform of a solid sphere. When normalized to a value of 1, it is given by where R is the radius of the sphere. The simulation of the Coulomb explosion provides the position of each layer as a function of time. The scattering factor for a layer with inner radius R 0 and outer radius R 1 is . 23 The scattering factor for the 1D spherical sample at time t is then written as ; , , ) . The mutual intensity matrix is then constructed from equation (2) by integrating over time.
A 7.5 nm radius particle with the composition of GroEL (pdb code: 1GRL) was simulated with a 10 keV, 40 fs x-ray pulse focused to 100×100 nm 2 spot size for two different pulse intensities. At the weaker intensity,5 10 19 W cm −2 , the onset of the explosion occurs toward the end of the 40 fs pulse as shown in figure 3. The weight for the second mode is a factor of 0.005 smaller than the first mode, which suggests a single mode approximation is valid. At high resolution, however, higher order modes have an increasing significance. Since the diffraction pattern has a dynamic range spanning several orders of magnitude, higher order modes can dominate the first mode at high resolution even when the weights for higher order modes are small. This can be visualized by the ratio, The plot of the R in figure 3 indicates that for5 10 19 W cm −2 around 5% of the total scattering intensity is due to higher order modes at resolutions beyond 0.15 nm −1 (i.e. beyond the first speckle). When the intensity is increased to5 10 20 W cm −2 , which is close to the maximum intensity available at LCLS [22], the explosion starts at around 10 fs. In this case, the ratio of the second modal weight to the first modal weight is 0.005, the same as for the weaker intensity. However, R is as high as 10% 20% everywhere except at very low q values. The realspace representation of the eigenvectors shows that the dominant mode resembles the original sphere with a reduced density at the surface for layers that move during the explosion. The second and third modes contain information about layers that are accelerated together during the explosion, producing correlated scattered intensity. This seems to indicate that the dominant mode provides information about the original structure and that higher order modes provide further information about the explosion.
For both intensities studied here, the first three modes account for almost all the scattered intensity at all resolutions. This is comparable in complexity to successful implementations of partial coherence CDI with a small number of dominant modes, suggesting that imaging an exploding particle may be possible using a similar scheme. The difficulty is that there will not be a priori information about the explosion to help constrain the reconstruction of higher order modes. Secondly, the 1D spherical modal only accounts for an average radial expansion. Further modes will likely be required for correlated motion in a 3D structure. The background due to uncorrelated motion, as described in section 4 would also be increased by the explosion. Nevertheless, it is interesting that the average radial expansion can be represented by only two or three modes.

Coherence in a spherical model of heterogeneity
From the early experiments XFEL single particle imaging heterogeneous samples have been identified a key issue. It can arise both as a byproduct of sample delivery or it can be intrinsic to the particle itself. Heterogeneity is problematic for 3D imaging experiments which require many copies of a particle to share a common structure. Data from a heterogeneous sample exhibits reduced contrast when averaged, which can be analyzed with methods of partial coherence. The resulting modes can be interpreted as dominant structural correlations in the heterogeneous sample. To illustrate such an analysis of heterogeneity, we have studied sphere diffraction data with variable radius. The mutual intensity function is constructed using the scattering factor for a sphere averaged over a distribution of sphere radii, as follows where P R ( ) is the probability of measuring a sphere of radius R. The distribution P R ( ) was chosen to be Gaussian with a standard deviation denoted by σ. Figure 4 shows the first three modes in q-space and real space for a mean sphere radius of 50 nm and standard deviations of s = 1 nm and s = 5 nm. In both cases, the first mode is the average sphere density, which tapers toward the surface in the region affected by heterogeneity. The higher order modes contain information about the surface. The weight of the second mode is 0.05% of the first mode for s = 1 nm and 1% for s = 5 nm, which indicates that in both cases the first mode dominates. As with the Coulomb explosion, the full story is only clear with a study of resolution. Figure 4 shows a plot of the number modes contributing more than 3% of the intensity of the first mode. It shows that the single mode approximation is good until the half-period resolution is equal to σ, at which point the number of modes required increases rapidly. This is consistent with the random heterogeneity model which contains no structural correlations. Such a rapid increase in modes suggests that using a partially coherent CDI scheme would provide little reward for effort. Heterogeneity effectively imposes a resolution limit.
The result for the heterogeneous sphere model helps to put the analysis of the spherical radial expansion in context. The radial expansion is not equivalent to a ensemble of random structures; as layers are sequentially accelerated away from the surface there are clear structural correlations captured by higher order modes.

Conclusion
Induced sample dynamics is an inevitable part of XFEL imaging and crystallography. Partial coherence analysis provides a way to understand the detailed implications of sample dynamics for structural analysis. The existing studies of electronic dynamics have provided the foundations for applying analysis in the context of XFELs and here we have identified extensions to crystallography, Coulomb explosion and heterogeneity.
In crystallography at the highest available XFEL pulse intensities, coherence effects impact the relative contributions of different atomic species. This is an important consideration for the development of anomalous dispersion and isomorphic replacement methods for XFEL crystallography, which use heavy atoms.
The impact on coherence of the average radial explosion of a particle, as represented by a 1D spherical model, was shown not to exceed more than two or three modes for a pulse intensity typical of current XFEL sources. As layers are sequentially removed by Coulomb forces, there remains a high degree of correlation in the model. Although the explosion dynamics of a 3D structure will be undoubtedly more complex than our 1D model, it would be worth exploring further the imaging of an exploding particle using partially coherent CDI schemes.
Finally, we have provided an example of coherence analysis for a dataset of spheres with variable radii. The scattered intensity can be described with a coherent scattering model, representing the average density, up to a resolution cut-off determined by the length scale of heterogeneity. In this case, the average diffraction does not retain information the fluctuations of the sphere radii, because these radii are uncorrelated.
The high-intensity, destructive pulses produced at highrepetition rates by XFEL sources deliver data that is both intrinsically statistical and impacted by radiation damage. Coherence analysis is both a powerful tool for understanding the diffraction of these pulses and for extending the capabilities of imaging techniques to extract more from experiment. There is potentially much to be gained from the further development of coherence analysis for XFEL science.