Strains and growth conditions
E. coli NEB 10-beta (NEB10β) was used for all plasmid cloning and protein production. For all cloning, E. coli strains were cultured in Terrific Broth (TB) containing 24 g/L yeast extract, 20 g/L tryptone, 0.4% v/v glycerol, 17 mM KH2PO4, and 72 mM K2HPO4 at 37 °C with appropriate antibiotics (50 µg/mL kanamycin). M9 glucose medium with tryptone supplement (2% w/v glucose, 1× M9 Salts, 75 mM MOPS pH 7.4, 12 g/L tryptone, 5 mM sodium citrate, 2 mM MgSO4 7H2O, 100 µM FeSO4 7H2O, 100 µM CaCl2 2H2O, 3 µM thiamine, 1× micronutrients [40 µM ZnSO4 7H2O, 20 µM CuSO4 5H2O, 10 µM MnCl2 4H2O, 4 µM H3BO3, 0.4 µM (NH4)6Mo7O24 4H2O, and 0.3 µM CoCl2 6H2O]) was used for protein production in bioreactors.
Chemicals and reagents
Unless otherwise noted, all chemicals and reagents were obtained from MilliporeSigma. Plasmid purification and gel extraction kits were purchased from iNtRON Biotechnology. FastDigest restriction enzymes and T4 DNA ligase were purchased from Thermo Fisher Scientific and used for all digestions and ligations following manufacturer protocols.
Construction and expression optimization of titin monomer and polymerization cassettes
The amino acid sequence of rabbit soleus titin domains I67-70 was obtained from a recent publication18, and the coding sequence was computationally optimized for E. coli expression using DNA 2.0 (ATUM) (Supplementary Table 4)30. The resulting optimized sequence was synthesized as a gBlock fragment by Integrated DNA Technologies. The sequence was then inserted between the KpnI and Kpn2I restriction sites of modified BglBricks63 vectors containing gp41-1C and gp41-1N SIs under the control of a PBAD promoter or a PLacO1 promoter, yielding plasmids p-1-4XT-1B and p-1-4XT-1L, respectively (Supplementary Table 5). Additionally, the optimized titin sequence was inserted between the KpnI and Kpn2I restriction sites of a modified BglBricks vector containing no SIs under control of a PLacOI promoter, yielding plasmid p-4XT (Supplementary Table 5). To construct the 8Ig titin plasmid, PCR was first used to amplify the optimized 4XT sequence from p-4XT, adding a Kpn2I restriction site to the 5′ end and maintaining a stop codon and a BamHI site at the 3′ end of the amplicon (Supplementary Table 6). This amplicon was then inserted downstream of the 4XT sequence in p-4XT via restriction digest and a two-part ligation, creating the p-8XT plasmid with the 4XT sequence duplicated (Supplementary Table 5). The 12Ig plasmid (p-12XT) was made in a similar fashion. First, the 4XT sequence was PCR-amplified twice, once with primers adding a Kpn2I site to the 5′ end and a SpeI site to the 3′ end; and another time with primers adding a NheI site to the 5′ end and maintaining a stop codon and a BamHI site at the 3′ end (Supplementary Table 6). The resulting PCR amplicons and p-4XT were digested with the corresponding restriction enzymes and ligated in a three-part reaction, yielding the p-12XT plasmid with a triplicated 4XT sequence (Supplementary Table 5).
Bioproduction in shake flask cultures
Overnight seed cultures of 50 mL TB medium were inoculated with single colonies carrying the desired construct (Supplementary Table 7). These seed cultures were then used to inoculate cultures of 500 mL TB in 2 L Erlenmeyer flasks at an initial OD600 of 0.08. Cultures were placed on reciprocal shakers at 350 rpm at 37 °C until OD600 reached 3.0, at which point the corresponding inducer was added (0.2% arabinose for p-1-4XT-1B and 1 mM IPTG for p-4XT, p-8XT, and p-12XT). Cultures were then continued at 37 °C for 20 h.
Bioproduction in fed-batch bioreactors
Both titin monomer and polymer were ultimately produced in 2 L fed-batch bioreactors (Bioflo120, Eppendorf). Transformants containing p-1-4XT-1L or p-4XT were cultured overnight in 50 mL TB medium at 37 °C on an orbital shaker. The overnight cultures were then used to inoculate an autoclaved 2 L Bioflo120 heat-blanketed bioreactor containing 1.5 L M9 glucose medium with tryptone supplement (see above). Antifoam 204 was added as needed to minimize foaming (approximately 0.01% v/v). Agitation and airflow were regulated to maintain approximately 70% dissolved oxygen (DO). After consumption of the initial 0.5% w/v glucose (as judged by ΔDO), a sterile substrate feed (20% w/v glucose, 48 g/L tryptone, and 10 g/L MgSO4 ∙ 7H2O) was initiated to maintain a linear growth rate. Reactors were induced at OD600 = 70 by addition of 1 mM IPTG, and the incubation temperature was reduced to 30 °C. Cultures were collected six hours after induction.
Protein purification
The 8 Ig, 12 Ig, and polymer titin were purified by resuspending cell pellets in urea lysis buffer (8 M urea, 300 mM NaCl, 10 mM imidazole, 20 mM KH2PO4, pH 7.4) at a ratio of 100 mL buffer to 50 g wet cell pellet weight. The solution was sonicated on ice using a QSonica Q700 sonicator (Qsonica) for 5 min (5 s on, 10 s off). Sonicated lysate was then pelleted by centrifugation at 25,000 × g for 30 min. Cleared supernatant was sonicated for an additional 5 min and then filtered through a 0.45 µm PES filter and applied to a series of six HisTrap HP 5 mL columns on an ÄKTA Pure Chromatography System (GE Healthcare Life Sciences) at a flow rate of 2 mL/min. Loaded columns were washed with two column volumes of lysis buffer, then washed by two column volumes of lysis buffer with 50 mM imidazole, and finally eluted by lysis buffer with 300 mM total imidazole. The ÄKTA chromatography system was controlled by and chromatogram data was acquired using the accompanying UNICORN software (Cytiva).
The titin monomer was purified by dissolving cell pellets in an aqueous lysis buffer (50 mM Tris, 50 mM NaCl, 1 mM PMSF, and 300 µg/mL lysozyme). After stirring for 30 min at 4 °C, 5 mM MgCl2 and 5 µg/mL DNaseI were added, and the mixture was sonicated with stirring on ice for 10 min (5 s on, 10 s off). After sonication, NaCl and imidazole were added to final concentrations of 300 mM and 10 mM, respectively. The mixture was centrifuged at 25,000 × g for 30 min at 4 °C, followed by 75,000 × g for 30 min at 4 °C. Cleared supernatant was then filtered and applied to a series of six HisTrap HP 5 mL columns at 2 mL/min. Loaded columns were washed with 2 column volumes of wash buffer (50 mM Tris, 300 mM NaCl, 10 mM imidazole), then washed with 2 column volumes wash buffer with 50 mM imidazole, and finally eluted with wash buffer with 300 mM imidazole.
After purification by affinity chromatography, the proteins were fully dialyzed to 5 mM ammonium bicarbonate at 4 °C using 10 kDa MWCO snakeskin dialysis tubing (Thermo Fisher Scientific).
SDS-PAGE
All SDS-PAGE gels were 1 mm thick, discontinuous with 3% stacking gel, and hand cast at the indicated percentages. Samples were prepared at 1 mg/mL total protein in Laemmli sample buffer (2% SDS, 10% glycerol, 60 mM Tris pH 6.8, 0.01% bromophenol blue, and 100 µM DTT). Gels were run on Mini-PROTEAN Tetra Cells (Bio-Rad) in 1× Tris-glycine SDS buffer (25 mM Tris base, 250 mM glycine, and 0.1% w/v SDS), until just before the dye front exited the gel. For MW estimation, we employed Precision Plus Dual Color Prestained Standards (Bio-Rad) and HiMark Pre-stained Standards (Thermo Fisher). Gels were stained in Coomassie Blue solution (50% v/v methanol, 10% v/v acetic acid, and 1 g/L Coomassie Brilliant Blue) for a minimum of one hour at room temperature with gentle agitation and destained in Coomassie Blue destain buffer (40% v/v methanol and 10% v/v acetic acid) for a minimum of one hour. Gels were imaged using the cSeries Capture Software on an Azure c600 Imager (Azure Biosystems). An unprocessed and uncropped image of the gel in Supplementary Fig. 2a can be found in the Source Data file.
Analytical SEC
Protein was concentrated to approximately 10 mg/mL based on absorbance at 280 nm. A Superose 6 Increase 10/300 column (GE Healthcare) was equilibrated with elution buffer (10 mM potassium phosphate, 150 mM NaCl, and pH 7.4), after which 100 µL of the sample were injected onto the column at 0.5 µL/min. The column was then eluted with 1 column volume of elution buffer, and the absorbance of the eluent was measured at 280 nm. Following the same procedure, 100 µL of protein standard mix (MilliporeSigma) and blue dextran (2000 kDa, MilliporeSigma) were separately passed through the column. A calibration curve was prepared by plotting the known MW of the standards against their retention volume (Vr) divided by the void volume (V0, blue dextran retention volume). An exponential curve was fit to the calibration data and used to calculate the MW of the titin polymer and monomer based on their measured retention volumes. Polymer number-average MW (Mn) was calculated as (1) ({M}_{n}=frac{mathop{sum}nolimits_{i}{M}_{i}{N}_{i}}{mathop{sum}nolimits_{i}{N}_{i}}), where (Mi) was taken as the calculated MW at a given data point on the polymer chromatogram (including only data from 1 kDa to 5 MDa) and Ni was taken as the measured absorbance at the corresponding data point. Weight-average MW (Mw) was calculated as (2) ({M}_{w}=frac{mathop{sum}nolimits_{i}{N}_{i}{M}_{i}^{2}}{mathop{sum}nolimits_{i}{N}_{i}{M}_{i}}).
Circular dichroism
CD spectra were acquired using the JASCO Spectrum Measurement software on a JASCO J-810 CD spectrometer equipped with a Lauda RM 6 refrigerated circulator and a JASCO PTC-423S peltier. Samples were diluted in 5 mM NaHCO3 and loaded into a 1 mm quartz cuvette (Hellma, Germany). The CD spectra were obtained at 20 °C, scanning in 1 nm steps from 190–260 nm with a 1 nm bandwidth, scanning speed of 100 nm/min, and 2 s response time. Multiple spectra were acquired, each the average of triplicate scans.
Scanning transmission electron microscopy
Samples (10 μL) were pipetted onto pure carbon, 400-mesh copper grids (Ted Pella, Inc.) that had been ozone-treated for 15 min using a Novascan PSD Series UV Ozone System (Novascan). After incubating for 5 min, grids were washed with 3 drops of ultrapure water and a 10 μl drop of 0.75% uranyl formate64. Grids were then stained by adding a 10 μl drop of 0.75% uranyl formate and incubating for 3 min. The filter paper was used to wick liquids away between each step, and the grids were allowed to air-dry before imaging. Samples were imaged on a scanning transmission electron microscope (STEM, JEM-2100F, JEOL, Japan) set at 200 kV using the GATAN DigitalMicrograph software. STEM images were simultaneously recorded from both a bright-field (BF) and a high-angle annular dark-field (HAADF) detector. Using ImageJ software (v. 1.52a), fibril cross-sectional diameters were measured approximately every 10 nm along the fibril axis, avoiding regions with ambiguous stain boundaries. A total of 376 diameter measurements were made.
Fiber spinning
Fiber spinning was performed by first dissolving lyophilized titin powder in hexafluorisopropanol (HFIP) to 20% w/v. This protein dope was loaded to a 100 µL Hamilton gastight syringe (Hamilton Robotics) fitted with a 23s gauge (116 µm inner diameter and 4.34 cm length) needle. The syringe was fitted to a Harvard Apparatus Pump 11 Elite syringe pump (Harvard Apparatus), and the dope was extruded into a water bath at 5 µL/min. Short segments (~5–10 cm) were then cut from these extruded fibers and carefully extended by hand in water at approximately 1 cm/s to 5× their original length. Extended fibers were removed from the bath and held under tension until visibly dry.
Light microscopy
Fiber diameters were measured using images acquired with a Zeiss Axio Observer ZI Inverted Microscope equipped with a phase-contrast 20× objective lens and the Axiovision LE software (Zeiss).
Scanning electron microscopy
Following tensile tests, titin fibers were mounted onto a sample holder using double-sided conductive tape (Electron Microscopy Sciences). The sample holder was sputter-coated with a 10 nm gold layer using a Leica EM ACE600 high vacuum sputter coater (Leica Microsystems). Fibers were imaged using a Nova NanoSEM 230 Field Emission Scanning Electron Microscope (Field Electron and Ion Company, FEI) at an accelerating voltage of 7–10 kV using the xT Microscope Control software (FEI).
Fourier transform infrared spectroscopy
For secondary structure determination, FT-IR spectra were acquired with a Thermo Nicolet 470 FT-IR spectrometer (Thermo Fisher Scientific) fitted with a Smart Performer ATR accessory with Ge crystal. Spectra were acquired from 1415–1780 cm−1 at 2 cm−1 resolution. A total of 254 scans were accumulated for each sample. All recorded spectra were analyzed using Fityk 0.9.865. Baselines were subtracted from all spectra using the built-in Fityk convex hull algorithm. The amide I band (1600–1700 cm−1) was deconvolved into a set of eleven Lorentzian peaks centered at 1610, 1618.5, 1624.5, 1632.5, 1642, 1651, 1659, 1666.5, 1678, 1690.5, and 1700 cm−1, corresponding to amide I shift characteristic of β-sheet, random coil, α-helix, or β-turn structures66,67,68. Peak areas were integrated, and component percentages were calculated as the component peak area over the sum of all peak areas. Percentages were averaged from measurements of three fibers for each condition (as-spun and post-spin drawn). To directly compare spectra, each individual spectrum was normalized to the highest measured absorbance. Normalized spectra were averaged (three spectra for each condition) and overlaid.
Polarized Raman spectromicroscopy
The method reported here is adapted from several previous studies of molecular alignment in spider silk fibers37,69,70. Titin fibers were carefully fixed to glass microscope slides with microscale markings to ensure that spectra were acquired at the same location before and after stage rotation. Raman spectra were acquired with a Renishaw RM1000 InVia Confocal Raman Spectrometer (Renishaw) coupled to a Leica DM LM microscope with a rotating stage (Leica Microsystems). Fibers were initially oriented along the x-axis (parallel to the laser polarization). Fibers were irradiated at a fixed point with the 514 nm line of an argon laser with polarization fixed along the x-axis and focused through a 50× objective (NA = 0.75). Spectra were recorded from 1100–1800 cm−1 with 1800 lines/mm grating. For each acquisition, a total of 10 spectra were accumulated, each for 10 s. The stage was then rotated to orient fibers along the y-axis with the same laser polarization, and spectra were acquired a second time at the same fixed point. No signs of thermal degradation were apparent, either visually or within recorded spectra. All recorded spectra were analyzed using Fityk 0.9.865. Baselines were subtracted from all spectra using the built-in Fityk automatic convex hull algorithm. For intensity ratio calculations, all spectra were normalized to the intensity of the 1450 cm−1 peak, which arises from CH2 bending and is insensitive to protein conformation69. For each fiber, the normalized intensity of the peak at 1670 cm−1 when oriented along the Y-axis was divided by the normalized intensity of the peak when oriented along the X-axis to give the intensity ratio (3) (I=frac{Y}{X}). This procedure was performed on a total of three separate fibers for each condition, and calculated intensity ratios were averaged. Spectra were also averaged and are presented with standard deviations for each point along with the spectra.
Wide-angle X-ray diffraction data collection
Synchrotron-based wide-angle X-ray diffraction (WAXD) analysis was performed on the BioCars 14-BM-C beamline at the Advanced Photon Source at Argonne National Laboratory, Argonne, IL. The wavelength of the X-ray beam was 0.886 Å, with fixed energy of 14 keV, and the beam size was 130 × 340 μm2 (horizontal × vertical). 2D diffraction images were recorded using a Pilatus3 S 2 M Detector and the samples-to-detector distance was 200 mm. CeO2 powder was used for the instrument calibration. For X-ray fiber diffraction measurements, the air background was measured first with no sample mounted on the sample stage. Then, bundles of 25 fibers, 1 mm in length and approximately 10 µm in diameter, were mounted across the opening of a rectangular paper frame. The assembly was loaded onto the sample stage with the fiber axis perpendicular to the X-ray beam, and the exposure time was 60 s to obtain a 2D diffraction image. The obtained diffraction intensities were subtracted by the air background intensity. Multiple images (≥3) were taken to improve the signal/noise ratio.
Wide-angle X-ray diffraction data analysis
To analyze the WAXD results, radial and azimuthal 1D profiles were sequentially obtained from the deconvolution of 2D diffraction images using the FIT2D software71. The deconvolution and fitting of 1D profiles were performed with the peak analyzer tool in the OriginPro 2016 software (OriginLab, Northampton, MA). The data were fitted with Gaussian functions using nonlinear least-squares fitting. We obtained 1D radial profiles of intensity versus scattering vector q (Å−1, radius within the 2D diffraction image) by integrating azimuthally over a sector typically 20–30 degrees wide along either the equator or meridian.
The WAXD analysis assumed an orthorhombic unit cell commonly applied to β-sheet crystallites in semi-crystalline fibers72,73,74. In particular, the 1D radial intensity profile along the equator includes two main equatorial Bragg reflections. Here the innermost equatorial peak is indexed as (200), corresponding to inter-sheet d-spacing along the unit cell a-axis and the outermost equatorial peak is indexed as (120), corresponding to inter-chain d-spacing along the unit cell b-axis (Supplementary Fig. 8)72,73,74. After deconvolution of the 1D profile, we obtained the peak center (PC), full width at half maximum (FWHM), and relative intensity (I) of crystalline peaks vs. amorphous components (Fig. 2h, i, Supplementary Table 2). The degree of crystallinity was estimated by dividing the intensities of crystalline peaks by the sum of intensities from crystalline peaks and amorphous components (Supplementary Tables 2, 3) (4) (% ({{{{{rm{Crystallinity}}}}}}=,frac{{I}_{{{{{{rm{Equatorial}}}}}}(200)}+{I}_{{{{{{rm{Equatorial}}}}}}(120)}}{{I}_{{{{{{rm{Equatorial}}}}}}(200)}+{I}_{{{{{{rm{Equatorial}}}}}}(120)}+{I}_{{{{{{rm{EquatorialAmorphous}}}}}}1}+{I}_{{{{{{rm{EquatorialAmorphous}}}}}}2}})). The center positions of the (200) and (120) crystalline peaks indicate a-axis inter-sheet d-spacing of 1.08 nm and b-axis inter-chain d-spacing of 0.46 nm, respectively (Supplementary Fig. 8, Supplementary Tables 2, 3). From the center position and FWHM of the (200) and (120) peaks, the Scherrer equation was used to determine the average crystallite size of 1.08 nm along the inter-sheet a-axis and 2.91 nm along the inter-chain b-axis, respectively (Supplementary Fig. 8, Supplementary Tables 2, 3). The Scherrer equation is expressed by (5) (D=frac{Klambda }{beta {{{{{rm{cos }}}}}}theta })75 where D is the mean size of the crystallite domains, K is a dimensionless shape factor (with a typical value of 0.9), λ is the X-ray wavelength (0.886 nm), β is the FWHM value in radians (conversion of the FWHM in our study to radians uses (6) (beta =2{{{{{rm{arcsin }}}}}}(frac{lambda times {{{{{rm{FWHM}}}}}}}{4{{{{{rm{pi }}}}}}})), λ is X-ray wavelength), and θ (°) is the Bragg angle). Conversion of the PC in our study to Bragg angle uses (7) (theta ={{{{{rm{arcsin }}}}}}(frac{lambda times {{{{{rm{PC}}}}}}}{4{{{{{rm{pi }}}}}}})times frac{360}{2{{{{{rm{pi }}}}}}}), where λ is the X-ray wavelength. Because the calculated average crystallite size along the inter-sheet a-axis is the same as the d-spacing, we suggest that β-crystals contain two β-sheets.
To estimate the degree of orientation of the crystallites along the fiber axis, we obtained two azimuthal 1D profiles from the 2D diffraction image (Fig. 2j, k). Figure 2j shows the plot of the diffraction intensity integration as a function of azimuthal angle within the radial range of the equatorial (120) peak, and Fig. 2k is the corresponding plot within the radial range of the equatorial (200) peak. After deconvolution of the 1D profiles, we applied the calculated FWHMs of the crystalline peaks and amorphous components to Herman’s orientation function (8) fcrystal = (3 < cos2φ > −1)/2 (Supplementary Tables 2, 3)39,40,41. Here φ is the angle between the c-axis and the fiber axis and < cos2φ > is obtained based on the equation (9) < cos2φ > = 1 − 0.8 < sin2(0.4 × FWMH(200)) > − 1.2 < sin2(0.4 × FWMH(120)) > . The parameter fcrystal is 0 for no preferred orientation and 1 if all crystallites are perfectly aligned39,40,41.
Fiber mechanical testing and cyclic loading measurements
Segments of post-drawn fibers (20 mm) were carefully laid exactly vertical across a 5 mm (vertical) × 15 mm (horizontal) opening cut into a 20 mm × 20 mm piece of cardstock and fixed with adhesive tape at both ends of the opening. Diameters of mounted fibers were then measured by light microscopy, averaging measurements at three points along the fiber axis. Mechanical properties were measured by axial pull tests on an MTS Criterion Model 41 universal test frame fitted with a 1 N load cell (MTS Systems Corporation). Cardstock holders were mounted between two opposing spring-loaded grips, and the supporting edges were carefully cut. Pull tests were conducted at a relative humidity of 45% and temperature of 22 °C, with a constant crosshead speed of 10 mm/min. Stress-strain curves were recorded by the MTS TestSuite TW Elite software using a 1 N load cell and a sampling rate of 50 Hz. Fiber breaks were recorded when a 90% drop from peak stress was detected. All mechanical properties were automatically calculated by the MTS TestSuite TW Elite software. Ultimate tensile strength was calculated as the maximum measured load over the initial fiber cross-sectional area (A = πr2), as determined from light microscopy diameter measurements. Modulus was calculated as the slope of a linear least-squares fit to the stress/strain data of the initial elastic region. Toughness was calculated as the area under the total stress/strain curve divided by the initial fiber volume (V = πr2h), as calculated from measured initial fiber diameters and set initial gage length of 5 mm. For each protein, a total of 14 fibers were measured in this manner.
For cyclic loading measurements, fibers were prepared and mounted as described above. Fibers were pulled at a rate of 10 mm/min to either: (1) a range of strains, starting near 0% and incrementally increasing until failure; or (2) a fixed, near-maximal elongation (12% for monomer and 30% for polymer fibers). They were then returned to 0% elongation, treated with 95% humidity air for 1 min, and pulled again. For each test, damping capacity was calculated as the ratio of the energy between the loading and unloading curves over the total energy under the loading curve. Damping energy was calculated as the energy difference between the loading and unloading curves divided by the initial fiber volume, as calculated from measured initial fiber diameters and set initial gage length of 5 mm.
Molecular dynamics simulation
The representative simulation volume of the fiber was constructed within the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS)76 with periodic boundary conditions to simulate bulk behavior. Visual Molecular Dynamics (VMD)77 was used for the visualization of our systems. Within the periodic box, we assembled an initial molecular structure based on our structural analyses of the real fiber. In detail, we first extracted the all-atomistic structure of the titin monomer (I67–I70, 4 Ig) from protein structure 3B43 in the Protein Data Bank18. The psf file was generated using the VMD “Automatic PSF Builder” toolkit. The four linked Ig domains form an oriented fibril. We aligned the 4Ig fibrils along the model y-axis and arranged 2 × 2 stacks of fibrils, with each 4 Ig fibril positioned antiparallel to its flanking 4 Ig fibril and offset by two Ig domains (Supplementary Fig. 12). Then we converted the pdb/psf files into a LAAMPS input file. Sodium ions were then added into the simulation box to neutralize the system. Considering that the titin fibers in our experiments were air-dried before pulling, we did not include explicit water molecules in our model. The CHARMM36 force field was used78, and the system included 23,312 atoms in total.
After minimizing the system energy for 10,000 time steps in LAMMPS, we first equilibrated the system in the NPT ensemble for 1 ns. Then at room temperature, we compressed the titin fiber in the x and z directions (perpendicular to the fiber axis) from 1 bar to 4000 bars and relaxed it back to 1 bar in 2 ns, followed by equilibration of 5 ns. The temperature and pressure damping parameters were chosen as 0.1 ps and 1 ps, respectively. The applied high pressure is to ensure that side surfaces of the fibrils have close contact with each other so that adjacent Ig-like domains can be paired to pack into the same crystalline domain as indicated by our WAXD results (Supplementary Fig. 8). The β-sheet content in the final model structure was measured to be 22%, in agreement with our FT-IR analysis. Further discussion of the model setup can be found in Supplementary Notes 3 and 4.
To carry out uniaxial tensile tests on the titin fiber, we conducted a non-equilibrium MD (NEMD) simulation in NPxPzT ensemble, with x and z directions at atmospheric pressure. Uniaxial strains were applied to the y-direction (fiber axis) with a constant strain rate of 1 × 108/s using the ‘fix deform’ command in LAMMPS and a time step of 1 fs. Simulations were carried out for 9 ns, with a final tensile strain of 90%. LAMMPS output the engineering strain and the stress, which consists of a kinetic energy term and the virial term (from interactions between atoms, such as pair, bond, angle, and dihedral contributions). All simulations were run three times.
To measure the intra- and inter-fibril bonded and non-bonded interaction energies, we input the static structures calculated from LAMMPS into the NAMD software79. Every structure was equilibrated for 100 ps and run for another 100 ps for energy calculation. We used the ‘NAMD Energy’ toolkit in VMD to output every intra- and an inter-fibril interaction term. The hydrogen bonds were calculated using the ‘Hydrogen Bonds’ toolkit in VMD, with only polar atoms (N, O, S) considered. The donor-acceptor distance and angle cutoffs were 3.8 Å and 30 deg, respectively. The salt bridges were calculated using the ‘Salt Bridges’ toolkit in VMD, with an O-N cutoff distance of 3.5 Å. The atomic stress was output from the ‘stress/atom’ command in LAMMPS, which also consists of the kinetic energy term and the virial term. During the simulation and measurement, the Lennard Jones cutoff was set as 12.0 Å, and the Particle Mesh Ewald method was used to account for electrostatic energy. Throughout all the simulations, the time step was 1 fs.
Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article.
