2.1.1. Simulated Systems

We considered a 23-residue long HAfps: wt GLFGAIAGFIEGGWQGMVDGWYG and E11A and and W14A mutants. The N-terminus was modeled as a charged amino group, and the C-terminus was amidated. For most cases, we simulated E11 in protonated (neutral) state; however, a peptide with its charged version, denoted as wt<sup>−</sup>, was considered as well for membrane-spanning configurations. Simulated systems included one peptide, 162 POPC molecules (81 per leaflet) and 9337 TIP3P [32] water molecules together with sodium and chloride ions necessary to construct a neutral system consisting of a membrane slab with ∼20 Å of 0.15 mol/L NaCl solution margins on both sides. Peptides and lipids were modeled with Amber99SB-ILDNP\* [33] and Amber Lipid14 [34] force fields, respectively. In addition, we considered wt HAfp simulations in transmembrane hairpin configuration using Charmm36 force field [35]. Starting geometries for surface bound and transmembrane peptides were taken from our previous runs [20]. Necessary mutations were introduced with Discovery Studio Visualiser (Biovia).

## 2.1.2. MD Simulations

MD simulations were carried out with Gromacs software [36]. They were constructed using periodic boundary conditions, with a time step of 1 fs (2 fs for CHARMM simulations), interatomic bonds constrained using LINCS algorithm [37], van der Waals interactions smoothly shifted to 0 at 10 Å (or cut off at 12 Å, with force-switch for CHARMM simulations), and electrostatic interactions calculated using particle mesh Ewald method [38] with 1.2 Å mesh spacing. Desired temperature and pressure of 1 bar were maintained by velocity-rescale thermostat and seimiisotropic ParinelloRahman barostat [39], respectively. Temperature replica exchange runs (tREMD) [40] were conducted using 24 or 40 replicas, at temperatures, *T*, ranging from 310 to 350 or 377 K, respectively (detailed list in the Supplementary Materials, Table S2). Exchanges were attempted every 1 ps. The diffusion of trajectories in temperature space was monitored to assure that each trajectory was able to reversibly sample the entire temperature spectrum. The summary of conducted runs is given in the Supplementary Materials (Table S3).

## 2.1.3. Kinetic Analysis

Opening and closing of HAfp structures was considered as a two state process. The assignment of peptide configurations visited during tREMD simulations to hairpin or boomerang states was based on root mean square deviation (RMSD) of backbone heavy atoms with respect to NMR hairpin structure (pdb 2kxa, model 1) [11], with dividing threshold of 2.5 Å that corresponded to a minimum in bimodal RMSD distributions obtained at *T* = 310 K. A set of temperature dependent kinetic equations was fitted to time evolution of open state fraction in tREMD trajectories and extrapolated to infinite time to give an estimate of equilibrium populations at *T* = 310 K, according to method introduced by van der Spoel and Seibert [41]. All calculations were carried out with the

use of g\_kinetics Gromacs module. Final fractions of hairpin structures were evaluated as an average of two asymptotic infinite time estimates based on two tREMD runs that started from fully closed or fully opened conformations.

#### 2.1.4. Free Energy Calculations

Potentials of mean force for peptide translocation along the axis perpendicular to membrane plane (*z* axis) were obtained using umbrella sampling simulations based on tREMD runs with 24 replicas. The peptides were restrained to hairpin geometry (reference NMR structure, pdb 2kxa, model 1) with harmonic potential, *Uh*, acting on pairs of C*α* atoms that were closer than 7 Å in the reference structure, with a force constant of 2.39 kcal/mol/Å2. Umbrella sampling was performed using biasing harmonic potentials with a force constant of 2.39 kcal/mol/Å2 that acted in the *z* direction between the center of mass of peptide 1–20 C*α* atoms (PCOM) and the membrane center (MCOM) defined based on the positions of three terminal carbon atoms of each lipid acyl chain, located within a cylinder of 30 Å radius and *z* axis passing through PCOM. PCOM relative positions, *zp* = *zPCOM* − *zMCOM*, were gathered for the replica run at T = 310 K. Window spacing along the reaction coordinate was 1 Å for *z* ∈ [−5, 16] Å range and 2 Å intervals for *z* ∈ [18, 34] Å.

To assess necessary equilibration length and the final time of production phase, *tend*, for each umbrella window, we determined the time, *tOK*, at which *zp* distributions gathered for *t* ∈ [*tOK*, 12 (*tend* + *tOK*)] and for *t* ∈ [ 12 (*tend* + *tOK*), *tend*] were similar according to Kolomogorov–Smirnov test with *p* > 0.1. Then, *zp* distributions gathered for *t* ∈ [*tOK*, *tend*], with the requirement that *tend* − *tOK* > 15 ns, were analyzed with weighted histogram analysis method [42] as implemented in Gromacs WHAM module, with a standard bootstraping error analysis.

The free energy cost of hairpin restraining with harmonic potential in membrane environment, *Uh*, was evaluated based on unrestrained runs using free energy perturbation formula [43]: Δ*G*0→*<sup>h</sup>* = <sup>−</sup>*kBT* exp(−*βUh*){*F*}, where *kB* is Boltzmann constant, *T* = 310 K. As a source of simulation frames, we used the last 500 ns of two tREMD runs (250 ns for E11A) that were started from closed and open structures (see above). A set {*F*} of 10,000 frames was randomly drawn from this pool such that to fulfil the proportion of open and closed structures as determined based on the kinetic analysis. The procedure was repeated 1000 times to obtain an average Δ*G*0→*<sup>h</sup>* and its error as a standard deviation. Analogous process in aqueous environment was split into simulation parts, in which the force constant of the restraining potential was gradually increased in steps 0, 0.001, 0.01, 0.1, and 1.0 to its full value used in *Uh*.

Free energy of unrestrained peptides as a function of PCOM position along the *z* axis, *zp*, defined as above was evaluated based on the probability distributions of *p*(*zp*) gathered during unrestrained tREMD simulations, according to *<sup>G</sup>*(*zp*) = −*kBT* ln *p*(*zp*) + *G*0 relation where *G*0 is an arbitrary constant. In the case of transmembrane configurations, an additional biasing potential was introduced: *Ub*(*zp*) = 12 *<sup>k</sup>*(*zp* − *<sup>z</sup>*0)<sup>2</sup> for *zp* > *z*0 and 0 otherwise, that prevented peptide from surfacing, with *z*0 = 6 Å for wt and E11A and *z*0 = 4.5 Å for W14A peptides. The resulting biased probability distribution, *<sup>p</sup>*(*zp*), was subsequently reweighted to obtain the unbiased *<sup>p</sup>*(*zp*), according to the following formula: *p*(*zp*) = *p*(*zp*) exp *βUb*(*zp*) *<sup>c</sup>*, with *c* being a normalization constant [44].

#### 2.1.5. Tryptophan Fluorescence Calculations

To estimate depth-dependent fluorescence quenching by brominated lipids based on our simulations, we adopted the model proposed by A. Ladokhin [45]. For a given average depth of lipid carbon-bound Br probe, *hm*, and its dispersion *σm* (both established based on pure POPC tREMD runs, see Table S1 for values), we calculated simulation averages of depth-dependent fluorescence profiles:

$$\frac{F\_m}{F\_0} = \sum\_{W \in \{14, 21\}} \left\langle \exp \left( -G(h\_W - h\_{\text{m}}, \sigma\_{\text{m}}, S) - G(h\_W + h\_{\text{m}}, \sigma\_{\text{m}}, S) \right) \right\rangle\_{\text{MD}} \tag{1}$$

where *Fm F*0 is the ratio of tryptophan fluorescence quenching by lipids specifically brominated at position *m*, to fluorescence without the quencher, *hW* is the position of tryptophan indole ring center along the *z* axis in simulation frames, and the two Gaussian terms, *<sup>G</sup>*(Δ*h*, *σ*, *S*) = *S σ* √<sup>2</sup>*π* exp (Δ*h*)<sup>2</sup> *σ*<sup>2</sup> , describe contributions from both membrane leaflets with *S* being the assumed quenching intensity [45] (see the Supplementary Materials, Figure S7 ). Given sets of *Fm*/*F*0 for Br probes at (4, <sup>5</sup>),(6, <sup>7</sup>),(9, <sup>10</sup>), or (11, 12) lipid acyl carbon atoms experimentally determined for each peptide, we used our simulation data to check what configuration of the respective peptide provided for the lowest root mean square error (RMSE) with respect to these values. To this end, we evaluated Equation (1) for the corresponding set of *hm* and *σm* values based on simulation frames representing the considered peptide configuration, and determined RMSE subject to *S* minimization. To assess the fraction *f* and 1 − *f* of transmembrane (TM) and surface (SURF) configurations, respectively, we assumed *Fm*/*F*0 = *f*(*Fm*/*F*0)*TM* + (1 − *f*)(*Fm*/*F*0)*SURF* and minimized the RMSE subject to *f* ∈ [0, 1] and *S*.

## 2.1.6. Membrane Perturbation

Lipid splays were defined as events when any of carbon atoms within lipid acyl chain was at least 2 Å further from membrane midplane than the phosphate atom of the same lipid. Lipids proximity to peptide was assessed based on the closest distance of their phosphate atoms to any heavy peptide atom. Lipids closer than 7 Å were considered as "close" and lipids further than 30 Å from peptides were considered for reference calculations. All results were block averaged, with block length of 50 ns, and the analysis was conducted for replicas simulated at 310 K.

Water membrane permeability, *P*, was estimated assuming inhomogeneous solubilitydiffusion mechanism [46], based on water density profile across the membrane, *ρ*(*z*), and position dependent water diffusion coefficient in *z* direction, *<sup>D</sup>*(*z*):

$$\frac{1}{P} = \int\_{z\_1}^{z\_2} \frac{\rho\_0}{\rho(z)D(z)} dz \tag{2}$$

with *ρ*0 denoting bulk water density, and integral extending between ±40 Å from membrane center (we note that *P* was rather insensitive to integration interval as long as it encompassed membrane interior; hence, we resorted to such simple choice). Water density profiles were evaluated from tREMD runs for replica at *T* = 310 K with gmx density tool. For the calculation of diffusion profiles tREMD trajectories were demuxed, and continuous trajectory fragments that remained at *T* = 310 K for at least 10 ps were used for analysis. The *z* axis was discretized into bins, *zi* of 2 Å width and if a water molecule was found within a given bin at time *t*, i.e., *z*(*t*) ∈ *zi*, its displacement within Δ*t* = 10 ps contributed to *<sup>D</sup>*(*zi*) = (*z*(*<sup>t</sup>* + Δ*t*) − *<sup>z</sup>*(*t*))2/2<sup>Δ</sup>*t*, where the averaging includes all such instances. To obtain final diffusion profiles for integration, *<sup>D</sup>*(*zi*) were interpolated with a Gaussian kernel, with *σD* = 1 Å dispersion. Water permeability obtained for pure POPC was (9 ± 2) × 10−<sup>3</sup> cm/s, with error based on three simulation blocks (30 ns each), which is in fair agreemen<sup>t</sup> with experimental value of (13.0 ± 0.4) × 10−<sup>3</sup> cm/s [47]. Water permeability in the presence of a peptide was calculated for the entire membrane patch and should be interpreted as corresponding to experimental conditions with 1:162 peptide to lipid ratio.

#### 2.1.7. Peptide Supervised Insertion

The procedure of peptide supervised insertion comprised 10 rounds of 50 simulations, 5 ns each. Starting structures included two random peptide placements and two selected from unconstrained simulations as already deeply inserted, all in hairpin conformation. After each round, the distance of W14 C*α* atom from membrane surface (center of mass of phosphate atoms within given leaflet) was evaluated within the 4 last ns of each simulation, and a frame with maximum insertion depth was selected as a seed for subsequent round. The probability of reaching particular position within *i* rounds was evaluated as *p* = ∏*i ni* 50 with *ni* being the number of rounds in which the maximum insertion depth was within 0.7 Å from the global maximum.
