The impact of radicals in cold atmospheric plasma on the structural modification of gap junction: a reactive molecular dynamics study

ABSTRACT Cold atmospheric plasmas are currently gaining increasing attention for cancer therapy. However, very limited studies regarding the interaction mechanisms between plasma species and tissues are available. We report the interaction of plasma produced species ( and ) with gap junction by employing reactive molecular dynamics simulations. Our results indicate and radicals can chemically react with N-terminal of gap junction resulting in its structural damage. There are two breaking mechanisms being identified: C-N peptide bonds and C-C bonds can be damaged by and radicals, respectively. Our findings could be particularly important for understanding the plasma-generated reactive species triggering bystander effect based on gap junction intercellular communication.


Introduction
Recently, the breakthrough in plasma has allowed it to be generated at room temperature and atmospheric pressure, so-called cold atmospheric plasma (CAP) [1][2][3]. CAP is a rapidly expanding emerging technology and has wide biomedical applications including sterilization, blood coagulation, wood healing, tooth bleaching, skin regeneration, and cancer therapy [4][5][6][7]. Main therapeutic effects are linked to presence of reactive oxygen species (ROS) and reactive nitrogen species (RNS) generated from the CAP [8][9][10]. ROS and RNS are known to regulate key biochemical pathways including physical and chemical changes within inter-and intra-cellular environments [11].
CAP is a surface or near-surface treatment modality; whatever the effects of the solvated ROS and RNS have in the surface layer of cells that are exposed to them, the effect on deeper layers of tissue must involve some cell-cell communication. On the other hand, the bystander effect is a biological phenomenon that refers to the process in which ionizing radiation causes death to cells directly exposed to the radiation whereas death of non-irradiated bystander cells would be provoked thereafter [12,13]. The bystander cells might be immediately adjacent (direct bystander effect) or certain distance away from the exposed cell (indirect bystander effect). In the latter case, the transmission of chemical signaling appears to be mediated by secreted soluble molecules. The candidates of signaling molecules can be nitric oxide gases, proteins, hormones or calcium ions [14].
Plasma generated ROS and RNS either enter the near-surface or act as ligands to membrane surface receptors and thereby initiate various cellular responses, including DNA damage, cell cycle arrest, and other redox-mediated stresses [9]. The stressed cells communicate with adjacent and distant cells through various processes and agonists, including cytokines. This cell-cell communication induces additional intracellular ROS creation or release which in turn causes subsequent stress responses in the local tissue [15]. On the other hand, cells experiencing DNA damage result in ROS levels to increase from mitochondrial sources, which initiates signaling cascades inducing gap junction (GJ) signaling to adjacent cells and the release of extracellular signaling compounds that the signal cell further away [12]. The signaling to other cells can take place from these bystander cells.
Although the underlying mechanism of bystander effect remains unclear, it has been demonstrated that direct irradiated cells might transmit the bystander effect to the neighboring cells through a process called gap junction intercellular communication [16]. Some studies have demonstrated that tumor cells are more sensitive than healthy cells to the bystander effect coming from irradiated cells, resulting in an advantage in tumor treatment [17]. As of plasma medicine, what the role gap junction plays in reducing or facilitating eradication of tumor cells upon the impact of plasma-generated ROS is not yet well understood. Therefore, the functional importance of gap junction needs to be investigated at molecular level. In this study, we first propose that plasmagenerated ROS induce bystander effect via gap junction, and we will confirm this hypothesis by using molecular dynamics simulations.
2 Structure of gap junction, molecular model and simulation method Gap junctions are specialized intercellular communication channels between neighboring cells with an apparent separation gap of 2-4 nm, which permit the exchange of various ions and small molecules through a regulated gate. A schematic picture of gap junction and connexins are shown in Figure 1. Each gap junction is composed of two apposed hemichannels termed as connexons contributed by each cell. Individual connexon is made up of six protein subunits called connexins (Cx) monomers formed from a family of 21 human proteins. The structure of connexin 26 (Cx26a typical connexin in humans) with molecular configuration can be obtained from protein data bank (accession no. 2Zw3) as shown in Figure 2. Connexin is a transmembrane protein with four transmembrane domains (TM1 to TM4). They are connected by two extracellular loops (EL-1 and EL-2) and one cytoplasmic loop (CL). Each connexin contains amino (NT) and carboxyl (CT) terminus in the cytoplasm.
Although multiple domains contribute to the biophysical functional roles of gap junction channels, certain of amino acids within N-terminal domains play important roles in unitary conductance, permeability/selectivity, and gating in response to transjunctional voltage [18]. Oshima et al [19]. reported the electron crystallographic structure of a human Cx26 mutant (M34A). They suggested that a bundle of NT can act as a 'plug' to close the channel by observing a prominent density within the pore of the channel on the 3D density map, which was supported by the subsequent structural investigation of a N-terminal deletion mutant of connexin 26 channel with decreased central pore density [20]. Recently, the crystal structure of the gap junction channel formed by human connexin 26 (Cx26) at 3.5 Å resolution was reported [21]. This structure does show that the N-terminal regions of the six connexins line the pore entrance and form a funnel structure stabilized by a circular network of hydrogen bonds [21,22]. Experimental studies also revealed that the main difference of three-dimensional  structure of connexin channels in a couple of states lie in a variety of conformations of the N-terminal domain [22]. Based on the structural and functional study, the pluggating model had been proposed and the general idea is that closure of gap junction channel is mediated by the blockage of N-terminal in the vestibule, and the head-tohead hemichannels can regulate their gating independently [22]. Therefore, the investigation of the interaction of plasma generated-ROS with gap junction (or more specifically, with N-terminal domains) would be essential for gaining a better insight in plasma-based cancer treatment.
The structure of N-terminal of connexin Cx26 is shown in Figure 3. Since Cx26 is a type of protein, it is a long-chain molecule consisting of a backbone made up of amino acids connected sequentially via a peptide bond. The interactions of RONS with such a protein result in protein and peptide oxidation, degradation, mutation, fragmentation and destruction, which all affect the overall stability and function of the gap junction.
ReaxFF [23] is a widely used, reactive force field for atomistic simulations of chemical reactions involving bond breaking and bond formation processes, offering fundamental atomic-scale insights into many biological reactions in a system. In ReaxFF, the bonded interactions are generated on-the-fly, based on distance-dependent bond-order functions. The potential energy can be expressed as which are bond, lone-pair, over coordination, under coordination, valence angle, double-bond valence angle penalty, valence angle conjugation, hydrogen bond, torsion, conjugation, van der Waals, and Coulomb energy, respectively. The detailed meaning of each term can be found in original references [23,24]. Bond order plays an essential role in ReaxFF, which distinguishes between contributions from σ bonds, π-bonds and double π bonds. Bond order needs to be adjusted to compensate for atomic over/under-coordination. The parameterized reactive force field covers more than half of the periodic table of elements and their compounds. Recently it has been extended to complex organic molecules such as DNA and glycine [25]. ReaxFF is optimized to obtain overall good accuracy similar to quantum mechanical calculations but at a much lower computing cost. As for modeling of plasma medicine, ReaxFF has been used in MD simulations to investigate the interaction of ROS with DNA [26], peptidoglycan and lipids [27], as well as liquid water [28]. In this study, we investigated the interaction of plasma species with N-terminal of connexin Cx26 by means of the ReaxFF reactive force field, especially using the C/H/O/N glycine/water parameters developed by Rahaman et al [25].
Our simulation setup is as follows. The N-terminal of Connexin Cx26 which consists of 20 amino acids was built in ArgusLab molecular modeling package [29] for 3D structure prediction. The geometry optimization was performed by initially choosing UFF (Universal force field) [30] and then PM3 after UFF. It should be stated that UFF is a typical classical or nonreactive force field, which was built in 1992. The potential energy is expressed as a sum of valence or bonded interactions and nonbonded interactions: which are bond, angle. torsion, inversion, van der Waals, and Coulomb energy, respectively. UFF (or other classical force fields) has fewer terms than ReaxFF since UFF cannot describe bond formation and charge transfer process. UFF covers all elements of the Periodic Table and it is suitable to describe the energy and dynamics of organic and inorganic compounds [31,32]. ArgusLab also offers higher accuracy of reoptimization using PM3 semiempirical level, which is used after UFF.
Then the optimized structure would be further used as the initial configuration for reactive MD simulations by LAMMPS package [33]. The dimension of the simulation box is 3.0nm× 3.0nm× 3.0nm where the N-terminal structure is placed in the center. The N and OH atoms in both ends are fixed to keep the structure from drifting. Periodic boundary conditions are applied in all three directions. In all simulations, the time step is set to be 0.1 fs. In the first stage before investigating the reaction between ROS and N-terminal, the structure of N-terminal is equilibrated at room temperature as follows. Firstly, N-terminal is relaxed under NPT ensemble (isotherml-isobaric) control for 100ps at 300 K and 0 GPa to obtain a stress-free structure, during which pressure and temperature are maintained in equilibrium with a heat bath and mechanical reservoir using a Nose-Hoover/Parinello-Rahman formalism [34][35][36][37][38]. Dynamic feedback between these reservoirs and the system allows the volume and the kinetic energy to fluctuate about the desired mean values as expected in an equilibrium system. The obtained structure is subsequently under NVT run for 50 ps. Finally, the structure is relaxed under NVE for 20 ps. In the second stage, 15 incident ROS are randomly placed more than 1 nm away from the N-terminal. The initial kinetic energy of incident ROS is set to be room temperature and their initial velocity directions are completely random. The detailed reaction mechanisms and pathways are carefully detected.

Results and discussion
The structure of the gap junction channel is formed by human connexin 26 (Cx26). Among those domains in connexin, above-mentioned N-terminal helix plays a critical role in controlling gating mechanism. The investigation of the interaction of plasma generated-ROS with gap junction (or more specifically, with N-terminal domains) would be essential for gaining a better insight in plasma-based cancer treatment.
In this study, we mainly focus on * OH and HO * 2 radicals. * OH radical is currently the most widely investigated species of ROS since * OH radicals as highly reactive species can also act as precursors for the formation of H 2 O 2 molecules. H 2 O 2 molecules do not directly interact with biomolecules but dissociatively react with each other, forming water molecules and HO * 2 radicals which can directly react with biomolecules. In this sense, H 2 O 2 molecules are believed to be more stable species in an aqueous environment and therefore are more likely to penetrate and interact with gap junctions.
Our results show that * OH and HO * 2 radicals do chemically react with N-terminal and can lead to the destruction of its structure. It was found that important bonds such as C-N peptide bonds can be broken by the impact of * OH and HO * 2 radicals. Mainly two breaking mechanism can be identified. It is observed that in most processes, the H-abstraction always takes place before subsequent reactions. The main difference should be noted that H atoms in the N-terminal are abstracted by * OH radicals whereas H atoms in HO * 2 radicals are abstracted by reactive site in N-terminal. The newly found C-N bond breaking mechanism is shown in Figures 4 and 5. As one HO * 2 radical approaches the OH group in Ser (see in Figure 4(a)), This OH group abstracts the H atom from the HO * 2 radicals, forming a water molecule and an oxygen molecule (see in Figure 4(b)), which is similar to the previous observation that the OH residue of the disaccharide in peptidoglycan can be abstracted upon impact of HO * 2 radicals [39]. After abstraction of OH group of Ser, a primary alkyl radical is produced in Ser, which is not stable. This finally leads to cleavage of C-N bond and subsequently creation of the double C-C bond and formation of a resonance-stabilized amide radical (see in Figures 4 (c) and Figure 5). The driving force for the C-N bond breaking reaction is that the resonance-stabilized amide radical is more stable than the primary alkyl radical. Such bond breaking reaction is termed as the homolytic cleavage. This process is very similar to the interaction of an * OH radical with Ala residue, during which * OH radical abstracts H atom in the Ala residue and results in cleavage of C-N bond and creation of the double C-C bond and formation of a resonance-stabilized amide radical in the same way [39].
It has been shown by MD simulations based on the self-consistent charge density functional tight-bonding (SCC-DFTB) method that three * OH radicals with two H-abstractions can lead to peptide bond cleavage [40]. This infrequent process can also be observed in reactive MD simulations as shown in Figures 6 and 7. After abstraction of the hydrogen atom of Gly by the first * OH radical, a secondary radical is produced in alpha-carbon, which is quite stable because the unpaired electron can be delocalized to the adjacent site [40]. Another *OH radical can occupy the same site to yield a hydroxylated alpha-carbon. The latter can react with other free * OH radical by Figure 5. Reaction mechanism of C-N bond breaking upon impact of an HO * 2 radical on OH group in Ser, with the creation of the double C-C bond and the formation of a resonance-stabilized amide radical. H-abstraction, forming a carbonyl group (see in Figure 6(a)) and finally yielding backbone oxidation and subsequent C-C bond cleavage (see in Figure 6(b)).
Recently, plasma undoubtedly has great potential application for cancer therapy.  [45]. HO * 2 is permeable across the mitochondria membrane directly and induce the leakage of high energy protons of proton-motive force. On the other hand, the ΔpH component of proton-motive force may provide a suitable chemical environment, and the H + from the redox proton pumps could be mainly localized on the surface of inner mitochondria membrane. Where, O • 2 − may consume the protons and conduct them across the membrane to the matrix side via a pK-dependent formation of a protonated species of HO * 2 (perhydroxyl radical) [46]. For endoplasmic reticulum (ER), the ER redox environment dictates the fate of entering proteins, and the level of redox signaling mediators modulates the level of ROS [47]. In the pathway of oxidative proteins, electrons transfer from protein disulfide isomerase (PDI) to molecular oxygen and endoplasmic reticulum oxidoreductin-1 (ERO-1) through a Flavin adenine dinucleotide (FAD)-dependent reaction. This electron transfer suggests that ERO1 has a strong association with protein load in the ER and can trigger ROS generation. Due to presence of transition metals (such as Fe 2+ or Cu 2+ ), H 2 O 2 can be converted to hydroxyl radicals ( * OH), which are highly reactive and can cause damage to lipids, proteins, DNA, and gap junction (GJ). Nitric oxide ( * NO) is a reactive radical produced from arginine by nitric oxide synthase (NOS) [48]. Nitric oxide has a very short half-life and can react with superoxide to form peroxynitrite (ONOO − ). ONOO − can modify the structure and function of proteins and is highly damaging to tumor cells. The HO * 2 radicals result in cleavage of C-N bond of Ser in Figure 4, and * OH radicals break C-C bond of Gly in Figure 6, resulting into GJ disfunction. ROS and RNS move into the bystander cell to induce cell death via by bystander effect, and then move in next cell, as shown in Figure 8

Conclusions
In order to elucidate the underlying reaction mechanisms, this reactive molecular dynamics simulation study demonstrates the possible interaction of * OH and HO * 2 radicals with N-terminal of connexon, which forms gap junction for intercellular communication. It is found that * OH and HO * 2 radicals do chemically react with N-terminal and can structurally damage it. Two breaking mechanisms are identified, involving C-N peptide bonds and C-C bonds by the impact of * OH and HO * 2 radicals respectively. It is confirmed that most processes are initiated by H-abstraction. The main difference between these two mechanisms is that H atoms in the N-terminal are abstracted by * OH radicals whereas H atoms in HO * 2 radicals are abstracted by reactive site in N-terminal.
The gap junction is highly exposed to ROS from plasma treatment. Due to the essential role gap junction plays in bystander effect, any structural change induced by ROS can affect its overall function and thus influence the effectiveness of cell-cell communication. We believe that this is crucial for plasma medicine based on ROS and we hope the research community could put more emphasis on the importance of gap junction intercellular communication under plasma treatment.
It has been recognized recently that more fundamental method such as self-consistent charge density functional tight binding (SCC-DFTB) method [49] especially DFTB3 method [50] is much more accurate and reliable than classical reactive MD simulations with regard to chemical reactions between plasma and biomolecules. In our future research, we will explore the same topic by using such more advanced technique to further validate current findings.

Disclosure statement
No potential conflict of interest was reported by the authors.