Introduction

Polymers are many-body physical objects; in order to describe their equilibrium state and dynamics, it is often required to translate chemical sequence information into free-energy landscapes in three-dimensional space. Rigorous theoretical descriptions can capture only special cases such as homopolymers or multiblock copolymers by following bottom–up approaches that start with interactions on the monomer level or considering the self-similarity of self-avoiding walks on the largest scales. The sequence space available by current polymer chemistry1,2,3 or in biopolymers exceeds the limits for closed physical descriptions and is not accessible for complete scans by molecular simulation techniques. A new paradigm of data-driven polymer science is increasingly encouraged by parallel sampling methods4,5 and the advances in machine learning (ML)6,7,8,9,10,11,12,13,14 and has the potential to explore yet undiscovered patterns in sequence–property relationships.

A prominent problem for sequence-controlled polymers is their transport through lipid membranes and biological barriers, which is linked to a wide field of potential biomedical and biotechnological applications. The translocation time of polymer chains through a narrow nano-pore on the scale of one monomer has been described for homopolymers15,16 by means of scaling relations and, later on, extended the theory to block copolymers17,18. As soon as local conformation entropy of the polymer comes into play by widening the pore to a finite diameter and length19,20, a general expression as a function of the sequence seems challenging in the moment for both charged and uncharged polymers. The absence of a closed analytic theory for sequence-controlled translocation meanwhile does not exclude technical applications of nano-pores for DNA sequencing21,22,23.

The picture is similar when considering the translocation of a polymer through a lipid membrane by direct penetration of the membrane’s core. Here polymer translocation can be considered as the diffusion of its center of mass along an effective free-energy landscape determined by the self-assembled membrane environment24,25,26. Translocation of homopolymers through bilayer membranes was recently described theoretically by means of propagators as the solution of Edwards equation27 in good agreement with coarse grained simulations28. Simulation results on random copolymers indicate that the main factors for copolymer translocation are their average hydrophobicity as well as their degree of adsorption29,30 at the membrane–solvent interfaces31, which shall be reflected in the main modes of their potential of mean force. Experimentally, the passive translocation of synthetic random copolymers into mammalian cells in the absence of cytotoxic effects was discovered32 and confirmed recently33,34. A rigorous theoretical description as a function of sequence, however, is missing to date. The lack of theory does meanwhile not exclude the recent progress in finding artificial cell-penetrating peptides and antimicrobial peptides by high-throughput screening35,36,37,38 that may even outperform evolutionary highly conserved Tat- or Penetratin-based sequences for biomedical application. Wimley et al. found that fine-tuned differences in short-block amphiphilic sequences have a significant effect on peptide translocation rates following rules that seems not obvious at the moment37. In turn, sequences leading to optimal points in their biomedical performance can be found in unexpected corners of sequence space that are potentially accessed by sequence–cargo co-evolution38. In this work, we make use of a massively parallel sampling of polymer conformations while scanning the full sequence space of a short copolymer and see that the intricacies of polymer sequence–property relationships can be subtle and unexpected already when considering relatively simple environments.

The laws of physics are, however, normally simple by means of requiring a relatively small number of parameters that can be extracted efficiently from high-dimensional data by ML methods. Artificial neural networks (NNs) have been successfully applied in the dimensionality reduction from chemical monomer composition of polymers to their material properties, such as glass transition temperature39,40,41,42,43, viscosity44, solvation free energies45, and electronic properties13,14 depending on the repeat units. When addressing long polymer chains, the training and test data are fundamentally limited to a fraction of sequences due to the exponential increase of sequence space. Any restriction or bias in the training data has yet undefined consequences for the NN’s projection into unexplored parts of sequence space. Efficient classification and optimization algorithms, such as based on artificial NNs, are in fact “black boxes” and their results therefore need better understanding and explanation. Theoretically, an NN can approximate any continuous mapping given that at least one hidden layer of neurons with sigmoid activation functions is contained46,47. Beyond that, the stacking of non-linear filters seems to mark a qualitative difference as compared to shallow ML algorithms such that they may develop internal representations of the input information that correspond to a hierarchy of abstraction levels. The distinguished generalization performance makes so-called deep neural nets particularly efficient when confronted with multiple tasks48 simultaneously, for instance, in finding quantitative structure–property or structure–activity relationships49,50,51,52. Recent advances in exploiting NNs for physical problems show that they can help to determine the essential order parameters necessary for predicting a mechanical state in future10 or classifying a magnetic phases8.

In this work, we have the luxury to access a complete sequence-to-property map available for training and testing NN algorithms thanks to the graphics processing unit (GPU)-accelerated4 sampling of random polymer configurations for a given sequence. GPU-accelerated Rosenbluth–Rosenbluth53 sampling of a copolymer in an external field modeling a lipid membrane allowed us to generate a significant number of configurations for all possible binary sequences for chain length up to N = 16. Based on this unbiased data, a NN is trained to predict mean first escape times of the polymer through the layer. By systematic selection of a training set, we tested the NN’s performance of projection into unseen parts of the complete sequence space.

The rest of the paper is structured as follows: In section “Results,” we introduce the sequence-complete sampling data set based on the Rosenbluth–Rosenbluth method for translocation time prediction. By comparison with free-energy estimates of self-avoiding walks near interfaces, we underline the physical meaning and richness of the results. We also analyze the performance of NN based on translocation time prediction for two different training schemes. In section “Discussion,” we summarize the results. In section “Methods,” we describe the polymer conformation sampling for estimating its translocation time through a membrane as well as the NN model applied.

Results

Rosenbluth–Rosenbluth sampling

Let us consider the inverse mean first escape time 1/τ as a measure for the frequency of translocation of a polymer through the membrane, which is presented in Fig. 1 as a function of the mean hydrophobic fraction along a backbone of N = 12 monomers. Results for all sequences are shown and grouped into point clouds centered at the corresponding ratios NT/N. The point clouds are shaped according to the number, nb, of blocks of H and T species along the sequence in a way that the points on the right hand side of a cloud represent a polymer with a larger number of blocks.

Fig. 1: Translocation frequency vs. fraction of hydrophobic monomers.
figure 1

Inverse mean first escape times (Eq. (7)) as a function of the fraction of hydrophobic monomers of a polymer of length N = 12. Results are shown for all sequences containing between two and nine T-type monomers. Results sharing the same number of T-type monomers are spread within windows of width 0.04 along the ordinate according to the number, nb, of H and T blocks within the sequence. The exact position along the ordinate is calculated as NT/N+0.04×[nb/N−1/2]. Results for seven sequences are highlighted by labels.

Full size image

The results in Fig. 1 confirm earlier predictions28 that a maximum of translocation frequency is found near a point of balanced hydrophobicity of the polymer as given by a balanced fraction of H and T units NT/N ~ 1/2, in case that the typical block size is in the order of the Kuhn segment of the polymer31.

In Fig. 2, we show the monomer sequences leading to the largest and lowest translocation frequency 1/τ as well as the results for triblock copolymers as a function of chain length for the balanced ratio NT/N = 1/2 of hydrophobic beads. For alternating sequences, the re-scaled translocation frequency (see Eq. (7)) remains in the same order of magnitude showing that the polymers are below the adsorption threshold for the given chain lengths. For polymers that are significantly localized at the membrane–solvent interface, one would expect that the desorption to be the rate-limiting process for translocation. Adsorption effect is clearly visible for diblock copolymers showing a nearly exponential decay of translocation frequency as a function of chain length. For diblock copolymers, we expect that the desorption of the hydrophobic block from the membrane is the most significant rate-limiting process, and consequently diblock sequences lead to minimal translocation frequencies. It is important to notice that, for hydrophilic blocks larger than the membrane width, the switch of a hydrophilic end from one solvent side to the opposing solvent does only require a limited number of hydrophilic beads to be in contact with the lipid core at the same time, whereas the escape of the hydrophobic block into the solvent requires all monomers of the block to be displaced into solvent environment. Dynamic barriers such as the steric hindrance of the polymer backbone by lipid tails, is, however, not included in the mean-field environment.

Fig. 2: Translocation frequency vs. chain length.
figure 2

Inverse mean first escape times as a function of the chain length for fractions of hydrophobic monomers of 1/2. We show results for sequences leading to maximal (minimal) inverse escape times.

Full size image

In Fig. 2, it becomes visible that the symmetry of the polymer sequence with respect to hydrophilic ends adds an important factor to the desorption probability, in particular when comparing results for triblock copolymers where the longest chains show a more than one decade larger translocation frequency as compared to diblocks. The difference can be understood qualitatively by estimating the adsorption free energy in the strong segregation limit as

$$\Delta {F}_{{\rm{ads}}}(N)=-c\epsilon {N}_{{\rm{T}}}+\Delta {F}_{{\rm{el}}}$$
(1)

where c is the average number of favored contacts a T monomer finds in the lipid environment (coordination number) and \(\Delta {F}_{{\rm{el}}}=-{k}_{{\rm{B}}}T\mathrm{ln}\,[{Z}_{{\rm{surf}}}/{Z}_{{\rm{free}}}]\) is an elastic contribution due to the reduction of the partition function from Zfree to Zsurf upon localization at the surface. The partition sum for a self-avoiding walk takes the form58,59,60

$$Z(N)\sim q{\mu }^{N}{N}^{\gamma -1}$$
(2)

where q is a non-universal amplitude that may depend on the particular form of short-range interactions and μ is the effective coordination number for the given random walk logic and lattice. The exponent γ depends on the topology of the polymer that is either in free solution or attached to a surface. One applies γ ≡ γ1 ≈ 0.67861,62,63 for strands having one end grafted, and γ ≡ γ11 ≈ −0.3962,63 for strands having both ends surface attached. The partition sum in free solution scales as \({Z}_{{\rm{free}}} \sim {\mu }^{N}{N}^{{\gamma }_{0}-1}\) with γ0 ≈ 1.156764,65,66. Since we further compare only ratios of partition sums for given total chain length, we assume that q- and μ-dependent contributions cancel up to a factor of the order unity.

The probability density to find a symmetric diblock copolymer in bulk solvent as compared to a state adsorbed at an interface as illustrated in Fig. 3 then reads

$${p}_{{\rm{di}}}(N)=\exp (\beta \Delta {F}_{{\rm{ads}}})={{\rm{e}}}^{-\beta c\epsilon N/2}\frac{{N}^{{\gamma }_{0}-1}}{{(N/2)}^{2({\gamma }_{1}-1)}}$$

Now, assuming that the desorption is the rate-limiting process, we write the estimate for the translocation frequency as

$${T}_{{\rm{di}}}={T}_{0}{p}_{{\rm{di}}}$$
(3)

In Fig. 2, we show the results for Eq. (3), where T0 = 0.123 and c = 19.6 have been adjusted for obtaining least-squared differences from the diblock Rosenbluth-Rosenbluth sampling (RS) results. The results confirm the dominance of the exponential factor resulting from pair interactions of the hydrophobic block.

Fig. 3: Free-energy profiles.
figure 3

Free-energy profiles for various polymer architectures (N = 12, NT = 6) and corresponding relevant states for estimating desorption probabilities.

Full size image

The ratio between partition sums for interface-adsorbed diblocks and triblocks allows to project from diblock to triblock predictions for translocation frequencies,

$${T}_{{\rm{tri}}}={2}^{2({\gamma }_{1}-1)}{\left(\frac{N}{2}\right)}^{1-{\gamma }_{11}}{T}_{{\rm{di}}}$$
(4)

which is plotted in Fig. 2 for comparison. The resulting up-shift catches up to the RS diblock results up to a factor corresponding to a remaining free-energy difference of 1.4kBT that is missing in Eq. (4). Note that we did not consider finite chain length effects in Eq. (2) in scope of this qualitative comparison.

With this discussion in mind, it is interesting to have a look back to Fig. 1 for understanding surprising features observed in the sequence maps of slightly hydrophilic polymers. By the example of a fraction of 4/12 of hydrophobic monomers, we demonstrate that the polymers comprising the shortest amphiphilic blocks (labelled by “(a)” and “(b)” in Fig. 1) are found in a middle range of translocation frequencies, while triblock copolymers similar as those discussed in Figs. 2 and 3 lead to the largest translocation frequencies. A comparison of the free-energy profiles shown as an inset in Fig. 1 underlines the interplay between surface adsorption and hydrophobic/hydrophilic balance that leads to the result. Polymers that contain short blocks (“(a)” and “(b)” in Fig. 1) are mainly subject to an effective free-energy barrier for insertion into the bilayers, which is the rate-limiting factor for translocation. The result reflects the fact that the polymer is effectively hydrophilic and shows negligible surface adsorption effects. Combining T-type monomers into a larger center block, however, allows for anchoring of the polymer at bilayer–solvent interfaces and thereby effectively reduces the rate-limiting repulsion from the membrane environment. On the other hand, for the diblock copolymers with NT/N = 4/12, adsorption at the bilayer–solvent interface turns over to dominate the free-energy profiles and leads to the largest escape times found for the given hydrophilic/hydrophobic ratio.

From the comparison between RS results for τ and previous work26,28,31 we therefore conclude that the dynamic interpretation of the sampling results is justified.

Machine-learned translocation times

The complete data set generated by GPU-accelerated RS sampling forms a powerful basis for benchmarking the ML-based search for sequences fulfilling given criteria. A network similar to Fig. 4 can be designed in order to predict sequence-determined properties of a polymer39,44,45. In this work, we stick to the example of logarithmic translocation times log(τ), and refer to a chain length of N = 14 monomers. The total number of sequences excluding the mirror-symmetric ones is S = 8256. The fraction of sequences within the training set we fix to ftrain ≡ Strain/S = 1/7. The ratio between training to the remaining test set results in 1:6. However, we follow two distinct schemes for the distribution of training sequences within the sequence space: In the uniform scheme, we define equidistant intervals of size 1/ftrain, along the τ-sorted sequences (id-space), and select the central sequences within each interval as the training set. In contrast, in the τ-window scheme we select every second sequence within a window S/2 < id ≤ S/2 + 2Strain, where “id” is a τ-sorted unique index in sequence space (section “Methods”). Note that thereby we select sequences within a narrow window in the upper half of translocation times.

Fig. 4: Architecture of the neural network.
figure 4

Neural network architecture for translocation time prediction of a polymer as a function of H/T sequence.

Full size image

In Figs. 5 and 6, we summarize the results of the training, and the performance of the resulting network with respect to the test set. In Fig. 5a, the development of the mean squared error (MSE) between NN- and RS-based \(\mathrm{log}\,(\tau )\) values for all test sequences (unseen) is presented. We note a reliable convergence of MSE values for both uniform and τ-window training sets toward a horizontal line indicating that training was stopped early enough for not running into over-training. In the case of the uniform training set, MSE results typically end up more than one order of magnitude lower as compared to the τ-window training set. The corresponding root mean squared deviation from the expected value typically reduces by a factor of \(\sqrt{2}8 \sim 5\). For the uniform training set, the root mean squared relative deviation from the RS-based \(\mathrm{log}\,(\tau )\) value points to a typical error of 1.0%, whereas for the τ-window training set we observe values of 7.3%.

Fig. 5: Performance of the neural network predictions.
figure 5

a Evolution of the mean squared error (MSE) of neural network-based \(\mathrm{log}\,(\tau )\) estimates for the test data set as a function of training epoch. b Relative error of the predicted value of τ according to Eq. (5) as a function of sequence index idtest in the test set. The result is averaged for 17 groups (bins) of sequences along the RS-τ sorted test set (idtest). Error bars denote a confidence interval of 90% by excluding the highest and lowest 5% of the data. The blue box labels the range of sequences of the training set in case of τ-window. c Percentage of correctly predicted faster or slower other sequences as a function of sequence index. Here we use the same binning and error bar definition as in b.

Full size image

In Fig. 5b, we show the corresponding mean relative error for the back-converted (not logarithmic) time τ according to

$$\frac{\Delta \tau }{\tau }=\exp [\Delta \mathrm{log}\,(\tau )]-1$$
(5)

where \(\Delta \mathrm{log}\,(\tau )\) is the absolute difference between the RS- and NN-based \(\mathrm{log}\,(\tau )\) values. For the uniform training set, the relative error scatters between −23% and +43% as found for the largest index idtest (largest τ), whereas for the fastest polymers 90% of sequences stay within an error of −6% to +15%. For this training set, equivalent to a random selection of sequences, such high accuracy of the network prediction is remarkable when seeing that the RS-based values of τ are spread by a maximum factor of τmax/τmin ~ 2 × 1011. For the τ-window training set, the relative error far away from the training window increases as compared to the uniform set. Nevertheless, as the maximum range of relative errors is found in the interval of −0.87 ≤ Δτ/τ ≤ 1.25 for the largest index idtest, we conclude that typically the prediction hits the right order of magnitude for τ despite the fact that we used only the narrow sequence window for training. It is interesting to note that the translocation times of the fastest sequences is typically predicted correctly by a factor of ~3 despite the large distance from the training window. Qualitatively, we expect that, when shifting the τ-window to smaller (larger) values of τ, the prediction accuracy for the smallest τ values will increase (decrease) while accuracy for the largest values of τ will decrease (increase), which is supported by preliminary data (not shown).

Absolute values are not always the main question for the modeled mapping; in some cases it is enough to obtain a decision statement upon the performance of two structures. When comparing two polymer sequences, for instance, we may ask which of those translocates faster. In Fig. 5c, we therefore show the performance of the trained network to give the right answer for this question as a function of sequence idtest. For this purpose, we calculated for each NN result τ1 for a given test sequence the fraction of all other test sequences leading to an NN output τ2 that holds the same relation τ1 > τ2 or τ1 < τ2 as the corresponding pair of RS sampling results. In case of uniform training, \(98.{7}_{-1.6}^{+1.1} \%\) of other sequences are correctly attributed as slower or faster (with a confidence of 90%), and for the τ-window training set \(96.{8}_{-4.5}^{+2.8} \%\) of pairs are correctly labeled. For the τ-window training set, the performance far away from the training window is reduced, in particular for sequences with a lower idtest index. However, the average fraction of correct decisions does not drop below 94.3% for the selected bin size.

By Fig. 5, we therefore demonstrated that a quantitative prediction of translocation times is possible by the applied ML model, and the accuracy depends crucially on the distribution of training sequences.

In Fig. 6, we outline more details of the training result by showing the predicted value of \(\mathrm{log}\,(\tau )\) for the whole test sets of uniform and τ-window in Figs. 6a and 6b, respectively. The monotony of the predicted data points for both training sets follows the base data line despite the scattering of the data as discussed for Fig. 5. In particular, for the τ-window training set, we emphasize that the order of translocation times is predicted correctly for the fastest sequences although the training set covers only a narrow window within the slower half of sequences.

Fig. 6: Unseen data predictions.
figure 6

Neural network prediction (dots) for the mean first escape time τ for unseen data (test set) are compared to RS-based results (gray line). Data are shown as a function of a unique identifier idtest for sequences in the test set that is sorted according to the RS-based result for τ. The ratio between training and test set sizes is 1:6. The number of hydrophobic units, NT, is shown as color-coded halos. In a, the uniform training set distributed homogeneously along the full RS τ-sorted sequence list. In b, we chose every second sequence within the blue labeled τ-windows range between idtest = 4128 and idtest = 5305. Training set sequences are skipped in this plot such that the slope is doubled as compared to the full sequence set within the labeled interval. The insets show details in the fields of lowest τ. The blue horizontal line in b inset labels the RS-based τ-value at the lower bound of the τ-window.

Full size image

Another interesting observation is the prediction of step-like features in translocation time (arrows in Fig. 6b) as function of idtest that are reproduced throughout the test set although located outside of the τ-window training range. Thus even the relatively simple network seems capable of finding a generic rule that links sequence and translocation time and thereby expresses the rather rich result based on Eqs. (7) and (8) without knowledge of conformation entropy nor the escape times. In view of the generalization performance observed for the τ-window training set, it therefore seems that the network developed an implicit internal representation approximating the mathematical rules linking copolymer sequence and translocation time.

Discussion

We apply a massively parallel sampling of the conformations of amphiphilic copolymers by means of self-avoiding random walks within a given density field representing a model for amphiphilic bilayer membranes. We estimated the free-energy profiles of the polymers composed of hydrophilic (H) and hydrophobic beads (T) with respect to distance from the membrane as a reaction coordinate. We calculated the mean first escape time τ as a measure for polymer translocation time through the model membrane all 2N binary sequences up to chain length N ≤ 16. Our results confirm that polymer translocation is controlled by a balance of the overall hydrophobicity of the polymer and is inhibited by adsorption at the bilayer–solvent interfaces26,27,28,31, which is consistent with the picture for small solutes67 and larger solid objects such as carbon nanotubes68.

Amphiphilic polymers at a balanced hydrophobicity show the smallest translocation times when the sequence exposes small repeating amphiphilic features, while longest waiting times are associated with a diblock structure of the whole chain. The different translocation rates between diblock and triblock copolymers as well as their chain-length dependence can be explained qualitatively when comparing adsorption-free energies at the bilayer–solvent interface involving surface-critical exponents. The relatively weak dependence of the translocation time of balanced hydrophobicity small-block alternating copolymers from chain length indicates that local amphiphilic features are only weakly interacting with the bilayer–solvent interfaces and the copolymer effectively resembles a homopolymer chain for which the membrane is energetically transparent. Chain-length dependence in this case is expected to increase when effective monomer association constants are stronger than in the present model. When considering slightly hydrophilic backbones, larger hydrophobic blocks start to become more prominent in sequences leading to smallest translocation times as they promote the association of the net-repulsive backbone with the hydrophobic membrane core.

The extensive database generated by RS sampling has been used to feed a multi-layer artificial NN with four hidden layers in order to explore the capability of so-called deep learning approaches for finding a general rule of how copolymer sequence translates into translocation times through biological barriers. The aim of this work is to test the meaningful interpretation of the “dirty work” of NNs provided by a complete data set of polymer sequences. We demonstrate that, even by using a low fraction 1/7 of uniformly selected training examples as compared to the total number 2N of binary sequences for N = 14, the NN achieves a root mean squared relative deviation in the order of 1% for the logarithmic mean first escape time \(\mathrm{log}\,(\tau ).\) In order to test the generalization performance of the network, we implemented a second training scheme, where training examples have been selected from a narrow window of sequences with respect to translocation times τ covering a factor of ≈30 between maximum and minimum translocation times contained in the training set. In this case, the network approximates the order of magnitude of the test data set covering a window being more than nine orders of magnitude wider. We conclude that the NN developed an internal representation of the mathematical rules linking sequence and translocation times, which involve a precise estimate of rate-limiting energy barriers. The network thereby encodes a complex interplay between polymer net hydrophobicity and sequence-dependent adsorption at the bilayer–solvent interfaces that to date can be treated in a theoretically closed form only for special cases as it involves the sequence-dependent polymer conformation entropy and solving the diffusion problem in inhomogeneous free-energy landscapes. Our results indicate a systematic decrease of prediction accuracy when moving into unexplored corners of sequence space and challenge future investigation on the relation between training data bias and prediction accuracy.

Methods

Rosenbluth–Rosenbluth sampling

We consider the diffusive transport of a polymer through a lipid membrane resembling a homogeneous oil slab as shown in Fig. 7. In particular, we are interested in mean first escape time of a polymer through the membrane as a function of length, N, sequence of hydrophilic head (H) and hydrophobic tail monomers (T). Coarse grained polymers are embedded into an external concentration field that represents bilayer membrane on a mean-field level composed of an hydrophilic region (H) and a hydrophobic core (T), as well as solvent (S). The hydrophobic core has a thickness of six lattice units.

Fig. 7: Simulation model.
figure 7

Illustration of a coarse grained polymer chains as used in our simulations within the grid occupancy of a laterally homogeneous membrane, and repulsive interactions between effectively two components (H,S) and (T).

Full size image

Monomers are represented as single-cell occupations on a simple cubic lattice, and bond vectors are taken from a set of 26 vectors with lengths of 1, \(\sqrt{2}\), and \(\sqrt{3}\) lattice units. Double occupancy of lattice sites is forbidden, and the monomers have excluded volume. This set of static rules corresponds to those of Shaffer’s Bond Fluctuation Model54.

Between hydrophilic sites (H and S), and hydrophobic sites (T), we implement short-range repulsive interactions. We write the internal energies of H and T monomers of the polymer as

$${U}_{{\rm{H}}}(\overrightarrow{r})=\epsilon {c}_{{\rm{T}}}(\overrightarrow{r});\ \ {U}_{{\rm{T}}}(\overrightarrow{r})=\epsilon ({c}_{{\rm{S}}}(\overrightarrow{r})+{c}_{{\rm{H}}}(\overrightarrow{r}))$$
(6)

where \({c}_{x}(\overrightarrow{r})\) are the number of lattice occupancies by species x on the 26 nearest neighbor sites55. In order to keep the model simple, we use only a single interaction parameter defined as ϵ = 0.1kBT with kB being Boltzmann’s constant and T the absolute temperature. For the enumeration of cx, the occupancy of the lattice by a given external concentration field (Fig. 7) is counted, and monomer–monomer contacts are taken into account in a way that contacts with the external field are screened by surrounding monomers. Thereby solvent-induced effects on polymer conformations are represented by the model.

For a given amphiphilic sequence, we aim to calculate the mean first escape time of a polymer between a repulsive boundary at z = −a and an absorbing boundary at z = +a, (Fig. 7)56,

$${\tau }=\frac{1}{D}{\int\nolimits_{-a}^{+a}}{\rm{d}}z{p}^{-1}(z){\int\nolimits_{-a}^{z}}{\rm{d}}z^{\prime} p(z^{\prime} )$$
(7)

where D is the diffusion constant of the polymer and p(z) is the probability distribution to find the center of mass of the polymer at a given distance, z, from the bilayer’s mid-plane. We define a = 22 lattice units, and D = 1(lattice unit)2, such that the dimensionless number of τ does not include the chain-length dependence of diffusion time.

The probability distribution p(z) is calculated by generating M polymer conformations \(\overrightarrow{R}=({\overrightarrow{r}}_{1},{\overrightarrow{r}}_{2},\ldots ,{\overrightarrow{r}}_{N})\) according to the RS scheme53. For each conformation \(\overrightarrow{R}\), the contact energy \(U(\overrightarrow{R})\) is calculated according to \(U(\overrightarrow{R})=\mathop{\sum }\nolimits_{i = 1}^{N}{U}_{X}({\overrightarrow{r}}_{i})\) in units of kBT according to Eq. (6) depending on the species X of the monomer X = H or X = T. The center of mass \(\bar{z}(\overrightarrow{R})=(1/N)\mathop{\sum }\nolimits_{i = 1}^{N}{\overrightarrow{r}}_{i}{\overrightarrow{e}}_{z}\) is evaluated with \({\overrightarrow{e}}_{z}\) being the lattice unit vector along the membrane’s normal direction. The distribution p(z) is then written as

$$\begin{array}{r}p(z)=\frac{1}{M}\mathop{\sum }\nolimits_{i,\bar{z}({\overrightarrow{R}}_{i})\simeq z}^{M}{W}_{i}{{\rm{e}}}^{-\beta U({\overrightarrow{R}}_{i})}\end{array}$$
(8)

where the condition below the sum illustrates that only those conformations contribute whose center of mass is found within a grid distance \((z-1/2)\;<\;\bar{z}\le (z+1/2)\) from z, and β ≡ 1/(kBT). In Eq. (8), Wi is the Rosenbluth weight of the ith conformation.

For a given sequence of H and T monomers in a polymer backbone, we calculate the mean first escape time according to Eq. (7) based on the generation of M = 1.5 × 107 RS-generated chains at uniformly distributed random positions within a periodic lattice of 64 × 64 × 64 lattice sites. The algorithm is implemented for GPUs4. In order to analyze how the mean first escape time depends on the amphiphilic sequence of the polymer, we perform the procedure for all 2N sequences for various degrees of polymerization N ≤ 16.

Artificial NN

We employ a fully connected NN involving \({\tan}{\rm{h}}\)-activation as sketched in Fig. 4. The network is composed by 2 hidden layers with 64 nodes each followed by 2 hidden layers with 32 nodes each. The input layer corresponds to a vector of values 0 and 1 representing the considered amphiphilic sequence of hydrophobic (0) and hydrophilic (1) monomers. The output layer consists of one neuron whose output is compared to the RS-based τ value for this sequence. The total network depth is n = 5, where only the output activation, \({\tan}{\rm{h}} [{\sum \nolimits_{i = 1}^{32}}({w_{n,i}}{h}_{n-1,i}+b)]\), includes a bias, b. Since absolute values of τ spread over several orders of magnitude, we perform the training with respect to its logarithm. The RS-based values of \(\mathrm{log}\,(\tau )\) are further linearly normalized and centralized into an interval [−0.9, 0.9] by defining \(I({\log}\,(\tau ))=1.8\times [({\log}\,(\tau )-{\log}\,({\tau }_{{\rm{min}}}))/({\log}\,({\tau }_{{\rm{max}}})-{\log}\,({\tau }_{{\rm{min}}}))-\frac{1}{2}]\) in order to be conveniently expressible by the \({\tan}{\rm{h}}\)-activation output. NN-based estimates for \(\mathrm{log}\,(\tau )\) are obtained by the back projection, I−1, of output neuron activations.

All weights are initialized with uniform random numbers in an interval [−0.3, 0.3]. The feed-forward (ff) back-propagation (bp)57 algorithm is employed for training. Error bp is performed after each ff cycle for a randomly selected sequence taken from the training set (stochastic gradient descent). The squared difference between the resulting activation of the output neuron and the RS-based \(I(\mathrm{log}\,(\tau ))\) value is used as the cost function for weight and bias adjustment. We set the initial training rate to η = 0.02, which gets reduced by a factor of (1/1.3) every 103 epochs in order to avoid frustration or early over-training effects. One epoch is defined as the average number of ff-bp cycles per sequence–τ pair. We set the total number of epochs to 104.

For each sequence, we define an unique integer identifier, 1 ≤ id ≤S, that is sorted according to the RS-based τ value. A lower id means a lower τ. The whole of S sequences is divided into a training set of size Strain and a test set of the size Stest = S − Strain. For the test set, we define a unique identifier idtest for each sequence that is the analog to id for the total sequence space. The index idtest labels sequences that are unseen by the network during training.