HIV-1 immature virion network and icosahedral capsids self-assembly with patchy spheres

The self-assembling process of a simple 3D molecular model that reproduces the hexagonal lattice that HIV-1 has in its immature stage is presented. The system is made of two overlapped repulsive spheres as a building block, both spheres are decorated with attractive patches simulating the main HIV-1 Gag interactions. If the overlapped spheres have the same diameter size, a flat lattice is assembled, while a curvature is induced in the resulting lattice when different diameters are used. Thus, the importance of this model lies on the fact that only by changing the overlapped spheres diameters it is possible to form icosahedral capsids according to the Kaspar and Klug model, known as T = 1, T = 3 and T = 4 without the use of any scaffold, such as the surface of a sphere. Only one molecular species is used to assemble hexagons and pentagons. WCA potential is used for repulsive particles and LJ for the attractive ones. Only two kinds of patches are needed to form a capsid. Two patches to form pairs, which form hexagons and pentagons, and the one that forms triplets, which connects the polygons together. These interactions are based on previous reports of the HIV-1 Gag protein. GRAPHICAL ABSTRACT


Introduction
Different types of hexagonal networks like Kagome and Graphene have been extensively studied, mainly due to the diverse applications that stem from the electric and magnetic properties of their structures.Graphene is obtained in the laboratory from graphite, but is also possible to computationally assemble these types of networks.One way to obtain these structures is starting with Janus particles [1] and colloidal multi blocks [2], to obtain the Kagome network [3].Besides, depending on the interaction sites positions and sizes, other network types may be obtained [4].
Fejer et al. [5] with a simple coarse-grained model involving anisotropic interactions can reproduce several morphologies, including spheroidal shells, tubular, helical and head to tail.Zhang et al. [6] used brownian dynamics to simulate large spheres composed of a smaller hard repulsive WCA spheres located on the large sphere surface, with different numbers and locations of LJ attractive patches also located on the large sphere surface.Several other geometries are also studied.Chains, sheets, rings, icosahedral, square pyramids, tetrahedron, and twisted and staircase structures are obtained.
Another significant hexagonal network is the one assembled by the HIV-1 virus.It is now known that HIV-1 Gag poly protein is formed by several individual proteins called dominions.From top to bottom these dominions are: Matrix (MA), Capsid (CA), Spacer Peptide 1 (SP1), Nucleo Capsid (NC), Spacer Peptide 2 (SP2) and p6 [7,8].Gag poly protein is responsible for the virus whole replication cycle.It is well known, since 1996, that the Gag poly protein MA dominion structures in a trimeric form [9]. Later on, Alfadhli et al. [10] showed these MA trimers forming a hexagonal arrangement over lipid membranes which has similar characteristics as the ones from the infected cells.Using cryo electronic tomography (CET) applied to immature virions it is possible to observe this hexameric arrangement presented by the Gag's CA dominion [11][12][13][14].Once the virus maturing process is terminated, a similar hexagonal network formed only by CA domain, is present in the capsid which protects the virus genetic material [15].During this encapsidation process hexamers have the capacity to become pentamers that help to completely seal the hexagonal structure.In this way, at the end, a conic capsid composed of 216 hexamers and 12 pentamers is assembled [15][16][17].
In addition to conical capsids, icosahedral ones are prevalent among viruses infecting humans, animals, plants or bacteria.There are capsids of different sizes; the smallest one form a regular icosahedral with 20 triangular phases, each composed of three proteins, with a pentamer in each of the 12 icosahedral vertices.A properly assembled capsid is built by 12 pentamers with a 60 proteins.Icosahedral geometry is kept for larger capsids formed by 12 pentamers connected with one or more hexamers.
Kaspar and Klug [19] developed an icosahedral capsid of different sizes classifying method.They assigned them a triangulation number T = 1, 3, 4, 7, 9, . .., where a larger number implies a bigger capsid.This classification is in accordance with the Goldberg polyhedra [20], characterised to be composed only by pentagons and hexagons, just like these types of capsids.This triangulation number provides information about the total number of proteins P = 60T and capsomers C = 10(T − 1) + 12; with 10(T − 1) hexamers and 12 pentamers.Triangulation number is defined as T(h, k) = h 2 + hk + k 2 , where h and k are positive integers that correspond to the number of hexagons to cross from one pentagon to the other within the same capsid.Traveling h hexagons in a straight line and k hexagons, to the right or left, to get to the next pentagon.
In the scope of numerical simulations, assembling of open networks in two dimensions, some of them with hexagonal geometries, is well studied.Using patchy particles, Kagome network has been assembled by Romano et al. [21] implemented the Kern-Frenkel potential [22] which imposes orientational restrictions between patches, and by Chapela et al. [23], without these restrictions.Edlund et al. [24,25] based on the uncertainty principle, used inverse methods to obtain isotropic interaction potentials that can assemble the Kagome network and several others.A variation of the Kagome called Deformed Kagome is reported by Machorro et al. [26] and the twisted Kagome by Salgado et al. [27] is assembled with patchy particles and a binary mixture of particles with isotropic potentials.A great variety of networks like Kagome, hexagonal honeycomb, square, rectangular, truncated square and truncated hexagonal have been assembled by Truskett group [28][29][30] using isotropic repulsive potentials obtained using inverse optimisation methods.
On the other hand, reports of these assembled networks in three dimensions are scarce, including the widely studied Kagome network.Fejer et al. [31] with a simple model of triblock Janus particles based on discoidal building blocks obtained Kagome structures in three dimensions.Romano et al. [21] assembled the Kagome allowing the spheres to rotate freely but restricting their displacement to the plane.One network that is easy to assemble in three dimensions is the hexagonal honeycomb.It can be assembled with particles having three attractive patches [32] or rigid trimers using the Reaction-Diffusion algorithm, based on the interaction distance of the extremes of the trimers without considering their orientations [33].
Edmundson [34] reported an early attempt to obtain a T = 4 capsid starting with a sphere having up to 60 charges in its surface.With respect to icosahedral capsids, self-assembly simulations are varied.Except for the T = 1 capsid, all other capsids are assembled using mixtures of various types, mainly pentamers and hexamers.All of them apply Kerner [35] rules or some variant of them.In a series of works, Rapaport [36][37][38] uses a cluster of several spheres to form a triangular block used to assemble various T = 1 capsids.With three different species of trapezoidal blocks, he assembles T = 1 and T = 3 capsids [39,40].Johnston et al. [41] assembled T = 1 capsid using pentagonal pyramids and T = 3 capsid with a mixture of pentagonal and hexagonal pyramids.Baschek et al. [42] also assembled T = 1 and T = 3 capsids with a repulsive sphere with three attractive patches around the equator to represent one protein, using 5 and 6 of these to form a pentamer and a hexamer.Capsids are then formed with a mixture of these two structures.Nguyen and Brooks [43] managed to correctly assemble up to capsid T = 13 starting with a mixture of pentagonal and hexagonal structures as constructive blocks.They interact with a series of specific rules proposed by Kerner [35].These rules change as the triangulation number changes.
An important effort to obtain the self-assembly of virus capsids, using seeds or mixtures of pentamers and hexamers, is included in the following discussion of reported works.VanWorkum et al. [44] used dipoles, quadrupoles and hexapoles to assemble several open networks and some nano tubes using a MC simulation method.The hexapole assembles an open hexagonal network like the ones reported in this manuscript.The formation of an icosahedral shell is also obtained.An icosahedral shell structure is formed.In addition, mixed nanotubes were formed using pentagonal and hexagonal 'seeds'.Gomez-Llorente et al. [45] used a method of free energy minimisation depending on the radius of curvature of the capsid.Mixtures of pentamers and hexamers interacting with a multipolar expansion of the LJ potential are used as models to obtain some capsid structures.Luque et al. [46] using a mixture of pentamers and hexamers represented by LJ particles, carried out MC calculations to obtain elongated capsids and T1, T3 and T7 icosahedral capsids.Fejer et al. [47] obtained shells from mixtures of pentagonal and hexagonal pyramids.
The next few references use different scaffolds to obtain closed structures with icosahedral symmetry with great success.Zandi et al. [48] used a LJ potential to perform MC simulations of two different kinds of capsomers, pentamers and hexamers, by varying well depth of each species.Simulations are carried out on the surface of a sphere where the molecules are allowed to change from pentamer to hexamer and vise versa.They obtain structures with T = 1, 3, 4 and 7, among others.Chen et al. [49] performing MC simulations studied a rigid cone-shaped model with 4 to 17 overlapped spheres, with various angles.These particles, interacting with LJ and square well potentials, are initially located on the surface of a hull sphere and allow to move on it.They report the formation of clusters of icosahedral symmetry corresponding to T = 1, 3, 7 and 13 virus structures.Pak et al. [50] used a coarse grain model with 161 beads per Gag, interacting through a long-range 12-6 LJ potential with a modified soft core are used to produce an open hexagonal lattice with the aid of a cell membrane and RNA.Various types of MD simulations with the LAMMPS [51] package are performed.They report the self-assembly of a hexagonal open network, formed on the cell membrane used as a scaffold.Guo et al. [52] Performed reaction-diffusion simulations on a coarsegrained model containing 6 sites over the surface of a sphere to assembly a hexagonal HIV open network.
In this work self-assembly of the immature HIV-1 hexagonal network formed by Gag poly proteins is presented.A simple model formed by patched overlapped repulsive spheres is used.The patches position and attractive interaction potentials are based on the principal actual interactions of the Gag CA domain: 2-Fold (2F), 3-Fold (3F), 6-Fold (6F) and IH [14,53,54].Figure 1 shows the selection process of the Gag poly protein main interaction sites to allocate the simple model patches used to self-assemble the network.One important result of this study is that changing one hard sphere diameter, it is possible to assemble icosahedral capsids.This is achieved without having an extra molecular species [42] or imposing interaction rules [43] other than the actual interactions of the HIV-1 Gag.Other significant and surprising characteristic of the simple models used is their capability to produce not only hexagons but also pentagons and other polygons as shown in Gutiérrez et al. [55].This characteristic allows the formation of closed capsids.
For the models used in this work, the 6-Fold interaction shown in panels (a,b) of Figure 1 is discarded.Even though, being one of the main interactions, it is (a) Immature HIV-1 hexagonal network construction with the CA domain using PDB 4USD [14] obtained from the Protein Data Bank (PDB) [18], which is a tool to explore the structure of 3D biological macromolecules.(b) Selection of only one hexamer signalling the sites where 2F, 3F and 6F interactions are located.The location of the intra hexameric interaction (IH) group is also indicated.Ip6 molecule is not shown in the 6F interaction, it is located at the centre of the hexagon and interacts with the six indicated sites.(c) Simplified representation of the hexamer substituting the MA domains by discs.(d) Final 2D model of the patchy molecule used to self-assemble the immature HIV-1 hexagonal network, excluding the 6F interaction.
between an external particle (Ip6) interacting with the six Gags forming the hexagon.The IH interactions, integrated with non-continuous regions at the base of CA domain.Lingappa et al. [53] considers that the multimerization by the IH's is also important.Models used here only three interactions are considered: 2F, 3F and IH.
The work is organised as follows.A Section on Models and Interaction Potentials is given next.A Simulation Method Section, where the methods used are described briefly, is followed by a Results Section where the HIV-1 Planar Hexagonal Network assembly, with and without 2F interaction, is described.The formation of icosahedral and non-icosahedral capsids is also included.A Conclusions Section closes the work.

Models and interaction potentials
In this section, description of models used, and definition of interaction potentials are given.

Models
Based on the work by Gutierrez et al. [55] in 2D, various models in 3D are developed in this work.Figure 2 shows the four models used to assemble the hexagonal network and the icosahedral capsids.All models are composed of a pair of overlapped repulsive spheres with diameters σ 1 = 1 and σ 2 ≤ σ 1 separated a distance d = 0.5 from their centres.These spheres are decorated with attractive patches around both their equators.The IH, 2F and 3F interactions of the CA domain that are mainly electrostatic are modelled with a LJ potential through two attractive patches of different sizes.Both IH and 2F interactions are represented with a small patch with diameter σ = 0.15 and are located at a distance δ = 0.425 from the sphere centre that contains them, in such a way that the patches are tangent from the inside of the sphere to which they belong.This size allows only the formation of pairs when they are united and interact with specific rules that will be given below.The remaining 3F interaction is represented by a larger patch of diameter σ = 0.45 located at a distance δ = 0.275 from the sphere centre of the containing sphere and remaining tangent to its interior.This patch only interacts with other patches of the same species, and they form triplets when they join.The position of the patches on the equator is as follows.The IH interaction is composed of two patches, with the same characteristics mentioned above, IH1 and IH2, located at an angle θ IH1 = 0 • and θ IH2 = 120 • , respectively.The remaining two patches 2F and 3F are situated at an angle θ 2F = 180 • and θ 3F = 270 • , respectively.Each of the models are described with more detail below.Spheres σ 1 and σ 2 in all models represent the two sub-domains that compose the CA domain: N-terminal domain (NTD) and the C-terminal domain (CTD).They are situated one on top of the other, as in the simple molecular models shown in Figure 2. In these sub-domains 3F interaction is located in the CA-NTD, while 2F and IH interactions are in the CA-CTD.However, to obtain the assemble of a HIV −1 network faster, the three interactions are in each one of the two central spheres.
Model M1, described in panel (a) of Figure 2, forms the hexagonal network without the 2F interaction.It is formed by two repulsive spheres of diameter σ 1 = σ 2 = 1 separated by a distance d = 0.5 between their centres.Both spheres include IH1 i , IH2 i and 3F i patches, with i = 1, 2 identifying the sphere to which they belong.Model M2 given in panel (b) of Figure 2 also has two repulsive spheres with diameters σ 1 = σ 2 = 1 separated at a distance d = 0.5 with the same patches as Model M1, but in this case, the 2F interaction on both spheres is included.Model M3, shown in panel (c) of Figure 2 has two repulsive spheres of different diameters σ 1 = 1 and σ 2 < σ 1 separated by a distance d = 0.5 between their centres.Different size diameters produce an angle φ between a straight-line tangent to both sphere's equators and a line tangent to the largest sphere's equator.This detail is shown in the frontal view of M2 in panel (c) of Figure 2.With an angle φ > 0, the molecular assembly shows a curvature towards the smaller sphere.Depending on this angle value, it is possible to obtain closed structures.The sphere with diameter σ 1 = 1 includes the IH1 1 , IH2 1 and 3F 1 patches.Again, in this model, the 2F interaction is discarded.For the smallest sphere σ 2 < σ 1 the same patches are used as in sphere 1, keeping the same positions on the equator, but their diameters are scaled with the diameter of their sphere σ 2 .Then the diameters of the patches are σ IH1 = σ IH2 = σ 2F = 0.15 σ 2 and σ 3F = 0.45 σ 2 .Similarly, their position δ is scaled to remains tangent inside to the sphere, remaining as δ IH1 = δ IH2 = δ 2F = 0.425σ 2 and δ 3F = 0.275σ 2 .Model 4 given on panel (d) of Figure 2 is identical to M3, with spheres of diameters σ 2 < σ 1 , but now the 2F interaction is added.

Interaction potentials
The spheres on a molecule are repulsive with others in different molecules, and only the patches on one molecule interact attractively with patches on other molecules.The interaction of two molecules is described by the following potential.
where U pp (r ij ) represents the attractive interaction between the patches given by: with np 1 and np 2 the number of patches of molecules 1 and 2 and with U LJ (r ij ) the LJ potential: with r ij and σ ij = (σ i + σ j )/2 the distance and mean diameter of the patches i and j and with ε ij = √ ε i ε j the geometric mean of the well depths of the potentials of patches i and j.
And U ss represents the repulsive interaction potential of the molecular spheres, given by: with ns 1 = ns 2 the number of spheres of a molecule and with U SLJ (r ij ) the WCA [56] potential displaced up with the cut off radius rc ij = 2 1/6 σ ij given by: Interaction between attractive patches is the same in all four models.Patches IH and 2F have values of the attractive potential parameters of well depths ε = 6, diameters σ IH1i = σ IH2i = σ 2Fi = 0.15 and cut off radius r c = 0.45.
While patch 3F has well depth ε = 6, diameter σ 3Fi = 0.45 with cut off radius r c = 1.Repulsive spheres of diameter σ 2 < σ 1 , the value ε = 6 is kept, but values of diameters and cut off radius are scaled with the diameter σ 2 of the sphere.Patch interaction rules for all the models are as follows.Interactions of sphere type 1 (S1) patches with sphere type 2 (S2) patches use a repulsive WCA [56] potential, given by Equation ( 5) with rc ij = 2 1/6 σ ij and σ ij the mean diameter of the two particles.There is only interaction between patches of the same type of sphere S1-S1 and S2-S2.Hence, when the patches belong to the same type of sphere the interactions are as follows: Patch 3F only interacts with itself, the same thing happens with patch 2F, interacting only with itself.IH patches belong to two different species IH1 (situated at 0 • ) and IH2 (situated at 120 • ) and they interact in a cross fashion: IH1 with IH1, while pairs IH1−IH1 and IH2−IH2 interact with a repulsive WCA [56] potential, given by Equation ( 5) with rc ij = 2 1/6 σ ij with σ ij the diameter of either IH1 or IH2.
For the vibrating models, a harmonic oscillator potential U B (r ij ), in units of ε, is used to represent the bonds between every pair of atoms.It is given by: where k = 500 is the potential constant in units of ε/σ , r0 in units of σ is the distance where the potential is zero, and r ij is the distance, in units of σ , between the bonded atoms.
The units for longitude, energy, and mass are, respectively: σ = 1, ε = 1 and m = 1.Hence, reduced units are defined as: distance r * = r/σ , energy U * = U/ε, density ρ * = ρσ nd and temperature Temp * = k b T/ε.With nd the number of dimensions; k b , Boltzman's constant; ρ = N/V, density with N the number of molecules and V the volume of the simulation box; Temp = 2Ke/k b ndN, temperature with Ke the kinetic energy.The value k b = 1 is taken, and reduced units are written without the asterisk for the rest of the text.

Simulation method
Molecular dynamics in the NVT ensemble is used in this work.All results are obtained using the simulation package HOOMD-Blue [57,58].In every simulated system, 512 molecules are used with a square initial configuration.A constant potential energy is used as a criterion for stopping the simulations.The step size is always the same and is fixed to dt = 0.001 in all the runs performed.HOOMD-Blue offers the option of using vibrating or rigid systems.Preliminary results with 256 molecules to assemble the HIV-1 planar network showed that the vibrating case produced a structure twice the size as the rigid case in the same simulation time.In the case of capsid production, capsid T = 1 assembled in nearly half the time for the vibrating case compared with the rigid one.Therefore, vibrating systems are chosen to perform all runs included in this work.To test the consistency of the results, a simulation with 1024 molecules is performed.The resulted configuration is very similar as the one done with 512, at the same temperature and density.

Results
Planar Hexagonal network production, Icosahedral and Non-Icosahedral capsids are the main results presented here.

HIV-1 planar hexagonal network
Two cases of assembly of planar hexagonal network are discussed: With and without 2-Fold interaction.

Without 2-fold interaction
Planar HIV-1 hexagonal network is successfully assembled with model M1.In the first two panels of Figure 3 a simulation result with 512 molecules at a temperature Temp = 0.66 and density ρ = 0.064 is shown.The initial configuration is a square array.In the upper panel the HIV-1 hexagonal network is clearly seen, with the hexagons painted in blue, to facilitate their identification, while molecules not forming part of any hexagon are depicted in red.The intermediate panel shows a lateral view of the assembled network, where it seems that several segments appear, but it is only one network composed by 55 hexagons traversing the simulation box due to the periodic boundary conditions.
Li et al. [32] assembled a good honeycomb hexagonal network in 3D using spherical particles with three attractive patches.The patches are symmetrically distributed within the sphere, and they only need to form pairs to assemble the network.On the other hand, the Kagome network assembled by Romano and Sciortino [21] with a one sphere model with two symmetrically distributed attractive patches forming only threesomes, is only accomplished by restricting the spheres to a plane.In contrast, the HIV-1 hexagonal network is formed by molecules where the attractive interacting patches are not symmetrically distributed within the spheres and the molecules need to form twosomes and threesomes.In this apparently more complicated model, a slightly elongated molecule, with a distinct interaction between upper and lower patches, stabilises the third dimension and propitiates a good assembly in 3D.
If a longer molecule is used based on the three principal Gag domains MA, CA, and NC composed of three tangent spheres each one having attractive patches, hexagons are equally formed.Several chunks of network are obtained without generating a single continuous structure.This option was not used in this work because it is a more complex model than the one, with two overlapped spheres, employed here.In the case of HIV-1 virus, even though the Gag poly protein is considerably long, due to its 6 domains, these problems do not arise.This is because the hexagonal assembly is produced when the Gags are attached to the lipid membrane of the infected host cell.When this happens, the Gags are aligned with each other over the surface of the lipid membrane, only diffusing on it.Other factor that can contribute to align the Gags is the RNA which is attached at the bottom side of the protein.Due to this fact, in a simulation where molecules are freely moving without any restriction, it is adequate to use a molecule with very well differentiated upper and lower parts.This is the case of all the models used here having two overlapped spheres.
Panel (c) in Figure 3 shows the hexagon counting along isochore simulation runs using model M1.Isochore density is ρ = 0.064 and the temperature varies from Temp = 0.58 to Temp = 0.7 with increments of Temp = 0.02.The plot shows clearly that only two runs, with Temp = 0.64 and Temp = 0.66, produce above 50 hexagons.In both temperatures, a maximum of 55 ± 1 hexagons are produced from 340 × 10 6 steps of the total 400 × 10 6 steps run.The next temperature Temp = 0.68 has a maximum hexagon counting of 46 by the end of the run.At this temperature, growth on the number of hexagons starts at 80 × 10 6 steps where it continuously increases, while the two best runs start at 10 × 10 6 steps.This growth not only starts later, but is increases slower than the best runs.It continues growing by the end of the run, but it is going to need a very long time to reach the best value.Simulation run at Temp = 0.7 ends with the production of only two separated hexagons.For the lower temperatures Temp = 0.58 to Temp = 0.62 hexagon formation has a good initial growth rhythm, but in any of the three cases reaches a value close to the 55 produced by the two best runs.In all three cases, hexagon counting reaches a maximum of 32, 48 and 47 hexagons for Temp = 0.58, Temp = 0.6 and Temp = 0.62, respectively.From the hexagon counting growth a small region of temperatures stands out, where the network formed reaches the largest possible size.Above this temperature region, growth slows and then reaches a limit where it stops, no hexagons are formed.Below this narrow optimum region, the network does not acquire its optimum size.

With 2-fold interaction
With model M2 including the 2-Fold interaction, hexagonal network is also assembled.Panel (a) of Figure 4 shows a hexagonal network with 57 hexagons using model M2.In this network, besides painting in light blue those molecules forming hexagons and in red or orange those not forming hexagons at the exterior, some molecules, not forming hexagons but delimiting holes at the interior of the network, are painted in dark blue.
Defects that appear preventing the formation of a continuous network are holes provoked by the absence of hexagons at the network's interior.This fact is in concordance with studies performed on the structure adopted by Gag proteins in immature virions [11][12][13].In these studies, only one interconnected hexagonal network is observed, having several defects caused by the lack of Gags and not for disordered Gags.These defects allow the network to bend, forming a spherical virion.In a hexagonal network, the inclusion of pentagons allows the network to bend and close to a spherical form.However, for the HIV-1 case Brigs et al. [12] conclude that the defects do not have a preference to adopt pentagonal forms, but they take diverse geometries, as is the case of the network in Figure 4.
Panel (b) in Figure 4 shows the hexagon counting for model M2 along an isochore at the same density ρ = 0.064 of the one calculated for model M1 shown in Figure 4. Due to the inclusion of the 2-Fold interaction, temperatures where there is a better network formation increase regarding those found for model M1.Temperatures on the isochore go from Temp = 0.76 to Temp = 0.9 with increments of Temp = 0.02.The results are like those of model M1.A broader range of temperatures where the network reaches a larger number of hexagons is found, going from Temp = 0.78 to a Temp = 0.86, with 54 to 59 hexagons.Above this optimal region, it is again observed a decay on the network growth rhythm until a point is reached where there is none.Below the optimal region, small groups of hexagons, when joined, do it incorrectly, without the possibility of rearranging to grow into a larger network.
Results for model M2 shown in Figure 4, inside the temperature region where the network stabilises forming a larger number of hexagons, for temperatures Temp = 0.78 and Temp = 0.8 the network is formed with defects.In temperatures, Temp = 0.82 to Temp = 0.86 formed networks are continuous and free of defects.A similar thing happens with model M1 shown in Figure 3, where in the temperature region of better network formation, the one formed at Temp = 0.64 has defects, while at Temp = 0.66 has none.The network growth process starts with the union of small hexagon groups formed at the beginning of the run, in the first few million configurations.Then the resulting network of this initial process continues adding hexagons to the periphery, composed of the free molecules available.It is in this initial growth process that defects are created.However, at higher temperatures, hexagons have a greater capacity to accommodate themselves, producing no defects in this initial network formation process.
Immature virions CET's [12,13,59] show that hexagons at the extremes of the network, including defects, have a lower internal hexagonal order than the other hexagons.From the processing and analysis of these data, Tan et al. [60] concluded that this is due because the extremes of the network of immature HIV-1 network are formed by incomplete hexagons.
Both networks formed with models M1 and M2 are surrounded by molecules painted in red joined to the periphery of the network.A similar thing happens with defects in the network's interior, where there are dark blue molecules, see Figure 4.In accordance with Tan et al. [60] these dark blue and red molecules can be identified as forming incomplete hexagons, lacking one or more molecules.Hence, the networks assembled in this work are in accordance with the observed immature virions.
In panel (a) in Figure 4 there are some molecules that are coloured in orange belonging to four incomplete hexagons with 5, 4, 3 and 2 molecules.They are located at the network extremes.Within the largest network defect, there are some incomplete hexagons marked.

Icosahedral capsids
To assemble icosahedral capsids models M3 and M4 are used.In these models, the two repulsive spheres have a diameter difference forming an angle φ between them, see panel (c) in Figure 2.This angle value dictates the degree of curvature of the hexagonal network.To obtain closed structures, it is necessary to have pentagons besides the hexagons needed to form the network.Models M1 and M2 may occasionally form a few pentagons at lower temperatures, when self-assembling a planar hexagonal network, despairing at higher temperatures.However, curvature imposed by angle φ stabilises the existence of these pentagons at higher temperatures.Molecular models M3 and M4 can form pentagons and hexagons, so the size of the capsid obtained depends entirely on the value of angle φ.
Figure 5 shows the capsids obtained with models M3 without the 2-Fold interaction.Hexagons are painted in a cyan colour, pentagons in black and molecules not belonging to any polygon are shown in red for the large sphere and green for the small one.T = 1, T = 3 and T = 4 are the capsids formed in this work.All of them show icosahedral symmetry.Following the Kaspar and Klug [19] classification method, the capsomers number of a capsid is given by C = 10(T − 1) + 12; with 10(T − 1) hexamers and 12 pentamers.Capsomer numbers of the three assembled capsids are: For T1 C 1 = 12 pentamers; for T3 C 3 = 20 hexamers and 12 pentamers, and for T4 C 4 = 30 hexamers + 12 pentamers.Angles required to assemble the capsids correctly are: φ 1 = 12.7 • , φ 3 = 7.4 • and φ 4 = 6.2 • , for T = 1, T = 3 and T = 4, respectively.Same Figure 5 shows a T = 7 half capsid with the hexamers and pentamers are positioned correctly, which means that, for an angle of φ 7 = 4.6 • , a complete T7 capsid could be formed.However, it was not possible in this work to close such a capsid with the correct number of capsomers.This might be since for triangulation numbers T(h, k) with k > 0 as T(2, 1) = 7, T(3, 1) = 13 or T(3, 2) = 19 the capsids are enantiomorphus, which means they exist in their left (levogyrous) and right (dextrogyrous) configurations.Namely, a capsid may be assembled in two ways, starting in a pentamer, turning left and advancing k capsomers until the other pentamer is reached or turning right and following a similar path.Not imposing a restriction in the simulation, to take one of these alternatives, two enantiomers may be combined avoiding the correct assembly of T7 capsid.
In order to find the correct φ T angles that are used to assemble the capsids shown in Figure 5, an initial angle φ i is selected to search around it.This initial angle is determined from manually assembling some contiguous capsomers based on a dodecahedron for T1 and a truncated icosahedral structure for T3.During this process, the molecule's smaller sphere diameter is adjusted to keep the molecule tangent to its neighbours.Once this is accomplished, the angle formed by the two spheres of a molecule is φ i .Many simulations are performed to assemble the capsids shown in Figure 5 to determine the correct angle φ T that produces the best assembled capsids.Once angles φ T are determined for model M3, the same values are used to assembly capsids with model M4.To determine T = 1 angle, a value of φ i = 12.7 • is used as a first guess, resulting in a properly assembled capsid.For T = 3 and T = 4 the initial angle is chosen to be in the range of angles calculated, starting with a second sphere diameter σ 2 = 0.86 up to σ 2 = 0.92 with an increment of 0.01.This corresponds to angles in the range of 5.1 • to 8 • .The resulting angles where a closed capsid is formed are for T = 3 φ 3 = 7.4 • and for T = 4 φ 4 = 6.2 • .In the case of T = 7, an initial angle of φ = 4 • is first used, looking for up to φ i = 5 • with intervals of 0.02 • .
Besides forming capsids with model M4, which includes the 2-Fold interaction, two other molecular models were used.These models are identical to models M3 and M4, with their 3-Fold patch having a well depth ε = 4 while the rest of the patches keep their original well depth value of ε = 6.This difference is motivated by the attraction forces measurement of 2-Fold, 3-Fold and IH interactions of HIV-1 Gag.These data are calculated with PDBePISA [61], a tool that calculates the Gibbs free energy of the inter-phases produced among the interaction sites of proteins.Results from these measurements indicate that the 3-Fold interaction is three times weaker than 2-Fold and IH interactions.Due to these results, the 3-Fold interaction ε value is changed.However, as models M3 and M4 2-Fold and IH well depth value is ε = 6, to maintain the twosome and threesome formation at the same temperature, the well depth value for 3-Fold interaction is fixed at ε = 4. Table 1 gives the number of molecules, densities, and temperature ranges used in the correct capsid assembly for the four models used.M3 dε (without 2-Fold) and M4 2Fdε (with 2-Fold) are the names that identify these two new models.From the data in Table 1 a few interesting facts emerge.In the sixth column, it is evident that as the capsid size increases, the range of temperatures where they are formed decreases.The last two columns of Table 1 indicate the temperature and number of time steps needed to assembly the fastest forming capsid.Several observations can be made on the results.One that is very clear indicates that as the capsid size increases the longer time is needed to assemble it.Furthermore, it is also notable that for all capsids formed models M3 dε and M4 2Fdε assemble them faster than their counterpart M3 and M4 models.Another remarkable fact is that for capsids T = 3 and T = 4 models M4 and M4 2Fdε , that include the 2-Fold patch, are formed from 1.2 to 2.1 times faster than models M3 and M3 dε , which do not include the 2-Fold patch.For the case of T1 the inverse case applies.Now versions without the 2-Fold patch, M3 and M3 dε , are 1.2 to 1.5 faster forming than those, M4 and M4 2Fdε , which include the 2-Fold patch.Lastly, for triangulation numbers T3 and T4, the fastest forming capsid temperature is very similar for all models.Again, T1 capsids behave differently, being model M3 temperature which is closest to other triangulation numbers.

Capsid stability test
Capsid T = 3 assembled with model M3 dε is used to test how stable assembled capsids are.From the creation temperature Temp = 0.61 the capsid is submitted to an increase or decrease in temperature.Taking the last configuration of the simulation run as staring point, six independent simulations of 50x10 6 time steps are performed at the following temperatures: 0.6, 0.5, 0.4, 0.3, 0.2 and 0.1.In all cases, the capsids remain assembled without any inconvenience.
Temperature increase is performed in two ways.First, applying temperature increments of size Temp = 0.02.The last configuration of the previous run is used as starting point and run for an extra 50x10 6 time steps.Second, thermal shocks are applied to the last configuration of the original simulation run and carried out 50x10 6 more time steps, as in the case temperature diminution.Temperature increases are applied in multiples of of size Temp = 0.01.For the thermal shock case, capsids kept their stable structure until Temp = 0.69, completely disintegrating at the next temperature.For the quasistatic increment, capsids maintain their stability until Temp = 0.73.

Non-icosahedral capsids
While looking for icosahedral capsids other types of closed structures without any holes are found.These constructs can also be called capsids, having 12 pentagons, as in the case of icosahedral capsids.The pentagons are not located at the icosahedral vertices but distributed in a distinct manner.The number of hexagons depends on the capsid size, which in turn depends on the molecular angle φ.
Figure 6 shows non-icosahedral capsids found.An interesting fact about these capsids, is that no matter the number of hexamers composing them, a certain symmetry in the arrangement of hexamers and pentamers is always kept.This is easily seen in panels (b) and (e) of Figure 6.Panel (b) shows an elongated capsid with a zigzagging ring of hexamers in the equator, while pentamers accumulate round the poles.In panel (e) a capsid in a tetrahedral shape is shown, presenting one hexamer in each of its four phases.These symmetrical arrangements of pentamers and hexamers are easier to identify in the first five capsids, but they all keep the same characteristic.
Panel (a) of Figure 6 shows T = 1 having only 8 pentamers.This capsid assembled two squares, one in each pole, to be able to completely close.On the one hand, this is made possible by the angle φ = 14 • used, which is greater than T = 1 capsid angle φ 1 = 12.7 • .This angle induces a larger curvature in the close structure formed.On the other hand, it is due to the model capacity to form various polygons, even when it is designed to form only hexagons.This capacity is based on the fact that the model central repulsive spheres can rotate around each other, maintaining the two IH forming pairs of patches at a distance equal to the minimum of their attractive potential, which is d min = 2 (1/6) σ IH = 0.1684.This rotating movement allows the angle to vary, if the two patches distance is always equal to the potential minimum, φ = ±11 • .This means that the angle to form pentagons (108 • ) and hexagons (120 • ), which are only 12 • degrees apart, are very close to φ, allowing the model to form pentagons and hexagons.

Conclusions
This work's main result is that it can reproduce large capsules from a simple model.Coarse-grained models with greater structural detail of the Gag CA domain of HIV-1 can assemble a 2D, hexagonal network [62].Even so, they do not obtain 3D closed structures, despite being able to form hexamers and pentamers [63].Furthermore, other simple models can assemble T = 3 capsids [42], but need to use three different molecular species.In contrast, the models presented here only require one molecular species as a building block.Therefore, changing the interaction rules to assemble a hexagonal lattice or an icosahedral capsid is unnecessary.
The models introduced in this work adequately form a 3D network, which leads to using them in future works, to try assembling other networks in 3D.Regarding the HIV-1 network, it is shown that the hexagonal structure can be formed without the need for 6F interaction, Ip6 molecules or other domains that can stabilise the virion.Another important fact is that the network can be assembled without the 2-Fold interaction.However, it is the only one that remains in the actual hexamers and pentamers of the capsid.This fact does not occurs with 3F and 6F interactions.
The appearance of pentagons on a hexagonal network allows to effectively close the structure.But the constructor block, in this case the molecules, must morphologically have the conditions to produce the curvature; this is to say, it must have a conical form.When the models presented here have an angle φ = 0 pentagons are produced but they are not stable.But if φ > 0, pentagons are adequately stabilised within the network.Measuring HIV-1 Gag dominions diameters, using the complete Gag reconstruction made by Sundquist et al. [7] Figure 1(c), a conical form appears with an angle φ Gag = 2.2.On the other hand, Pak et al. [50] using a coarse grain model with 161 beads per Gag finds that when the CA domain forms hexamers during the immature stage, the hexagon diameter obtained by its CA-NTD is larger than the one formed by its CA-CTD.Again, a conical shape is shown.They [50] also conclude that the CA domain induces a curvature on the hexagonal lattice, a fact that helps to bend the cell membrane that will form the HIV-1 virion.Thus, the models introduced here have the capacity to produce pentagons and hexagons with the same structural unit, along with the molecular model conical form, makes the formation of icosahedral capsids possible.Therefore, the only variable that is necessary to select is the angle φ, in order to produce one or other capsid.
The impossibility of the models used to form capsid T7 is perhaps due to the lack of control on which enantiomer must be assembled.Besides, in simulation works where enantiomers are chosen and polygons are used as capsids constructing blocks, the effectiveness of T 7 capsid formation is reduced below 10% [43].With the results presented here it is recognised that smaller capsids (T 4) can be assembled freely, but larger capsids require scaffolds [64] for a correct assembly.For the HIV-1 virus, these scaffolds may be their own genetic material and the lipid membrane of the infected cell [50].These elements may give the necessary support for assembling, firstly, the immature virion and subsequently capsids of a larger size.Many of the cited works [48][49][50]52] included, in their communications, one or both types of scaffolds, even for the self-assembly of T1 to T4 capsids.This is another theme for future work.

Figure 1 .
Figure 1.Attractive patch model development based on the HIV-1 Gag principal interactions that forms the immature network.(a)Immature HIV-1 hexagonal network construction with the CA domain using PDB 4USD[14] obtained from the Protein Data Bank (PDB)[18], which is a tool to explore the structure of 3D biological macromolecules.(b) Selection of only one hexamer signalling the sites where 2F, 3F and 6F interactions are located.The location of the intra hexameric interaction (IH) group is also indicated.Ip6 molecule is not shown in the 6F interaction, it is located at the centre of the hexagon and interacts with the six indicated sites.(c) Simplified representation of the hexamer substituting the MA domains by discs.(d) Final 2D model of the patchy molecule used to self-assemble the immature HIV-1 hexagonal network, excluding the 6F interaction.

Figure 2 .
Figure 2. Molecular models that form the immature stage of HIV-1 hexagonal network are shown on panels (a) and (b), and those that form icosahedral capsids on panels (c) and (d).Attractive patches represent the IH, 2F and 3F interaction of the CA domain of the HIV-1 Gag.The patches are located at the equator from each sphere at an angle 0 • and 120 • for IH1 i and IH2 i , respectively, and at angle 180 • and 270 • for 2F i and 3F i , respectively.The subscript i = 1, 2 distinguishes both spheres of the molecule.Models M1 and M3 do not include 2F interaction, while models M2 and M4 do.Models are shown with a certain opacity to allow the visualisation of the attractive patches tangent to the inside of their respective repulsive spheres of the molecule.

Figure 3 .
Figure 3. Hexagonal network of the HIV-1 immature virion assembled with model M1 without the 2-Fold interaction.Panel (a) shows the view from above of the network composed of 55 hexagons with parameters: ρ = 0.064, Temp = 0.66 and np = 512.Panel (b) gives the lateral view of the same network.Panel (c) shows the hexagon counting along the 7 runs at the same density ρ = 0.064 and particle number np = 512 with temperatures ranging from Temp = 0.58 to Temp = 0.70.

Figure 4 .
Figure 4. Hexagonal network of the HIV-1 immature virion assembled with model M2 including the 2-Fold interaction.Panel (a) shows the view from above of the network composed of 57 hexagons with parameters: ρ = 0.064, Temp = 0.8 and np = 512.Blue denotes molecules delimiting the inner network defects.Orange indicates some incomplete hexagons in the outside limits of the network.Black lines mark some incomplete hexagons in the network extremes and defects.Panel (b) shows the hexagon counting along 8 runs at the same density ρ = 0.064 and particle number np = 512 with temperatures ranging from Temp = 0.76 to Temp = 0.90.

Figure 5 .
Figure 5. Icosahedral capsids assembled with model M3.Density is ρ = 0.064 for all cases.(a) Capsid T = 1 formed only by 12 pentagons, with 60 molecules.(b) Capsid T = 3 is composed of 12 pentagons and 20 hexagons with 180 molecules.(c) Capsid T = 4 is composed by 12 pentagons and 30 hexagons with 240 molecules.These capsids are accompanied by a closed structure that has not the correct distribution of hexagons and pentagons.(d) Capsid T = 7 is not completely assembled, composed of 9 pentagons and 49 hexagons with 339 molecules.The required numbers are 12 pentagons and 60 hexagons with 420 molecules.Table I includes the number of molecule values, densities, and temperature ranges for the four models used.Capsids T = 3 and T = 4 show icosahedral symmetry.

Figure 6 .
Figure 6.Non-icosahedral capsids assembled with model M3.(a) smallest capsid found, formed by 8 pentamers and two squares, one in each pole, painted in red.(b) to (i) capsids of different sizes assembled with a diverse number of hexamers, always maintaining 12 pentamers.In all cases molecular angle φ is shown.