{'fraction_non_alphanumeric': 0.03129270673546228, 'fraction_numerical': 0.015032770882722075, 'mean_word_length': 5.049603509363928, 'pattern_counts': {'":': 1, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Can biophysical models give insight into the synaptic changes associated with addiction?

Effective treatments that prevent or reduce drug relapse vulnerability should be developed to relieve the high burden of drug addiction to society. This will only be possible by enhancing the understanding of the molecular mechanisms underlying the neurobiology of addiction. Recent experimental data have shown that dendritic spines, small protrusions from the dendrites that receive input from excitatory neurons, from spiny neurons in the nucleus accumbens exhibit morphological changes during drug exposure and withdrawal. Moreover, these changes relate to the characteristic drug-seeking behavior of addiction. However, due to the complexity of the dendritic spines, we do not yet fully understand the processes underlying their structural changes in response to different inputs. We propose that biophysical models can enhance the current understanding of these processes by incorporating different, and sometimes, discrepant experimental data to identify the shared underlying mechanisms and generate experimentally testable hypotheses. This review aims to give an up-to-date report on biophysical models of dendritic spines, focusing on those models that describe their shape changes, which are well-known to relate to learning and memory. Moreover, it examines how these models can enhance our understanding of the effect of the drugs and the synaptic changes during disease progression.


Introduction
The highest burden of drug addiction to all societies is the immense loss of human capital -from unfulfilled personal ambitions to loss of family structure. According to the National Institute on Drug Abuse (1 ), in 2019 more than 70,000 Americas died from a drug-involved overdose, including illicit drugs and prescription opioids. However, this number only represents a small percentage of drug-related deaths. According to a statistical study (2 ), in 2016 the drug-associated deaths was 2.2 times the number of drug-coded deaths. Besides these losses, the estimated medical cost in hospitals of substance use disorder in 2017 was $13.7 billion (3 ).
A continued medical aspiration for relieving the high burden of drug addiction in society is to develop and tailor effective treatment to prevent relapse. Assuming that drug addiction is a neuropsychiatric disease whose behavioral pathology consists in the propensity to relapse, even after periods of abstinence, then an effective treatment should prevent or reduce relapse vulnerability (4 ). Because current behavioral and pharmacological therapy only helps a small percentage of patients (4 ), we need to understand the molecular mechanisms underlying the neurobiology of addiction. In the past decade, studies have focused on the neurons in the nucleus accumbens (NAc), a brain region that mediates goal-directed behavior, and their structural changes after drug exposure (4 -7 ). These studies have noted that structural changes in neurons could be an important readout of addictive behavior. For example, when animals are exposed to an environment associated with previous cocaine consumption, it triggers cocaine-seeking behavior that is induced by structural changes in NAc neurons (7 ). Another study showed that structural alterations to these neurons are only present in animals that show a persistent increase in the psychomotor activating effects of cocaine after repeated exposure to it, i.e., psychomotor sensitization (8 ). However, there is no consensus on whether and how the structural changes in neurons derived from drug exposure drive addictive behavior (4 , 5 ) due to opposing experimental data which show that the structural changes in NAc neurons after cocaine exposure are not required for locomotor sensitization and suggest that, instead, these changes may represent a compensatory process (9 ). Complicating matters, studies have also shown that different drugs have different effects on these neurons (4 , 5 ). For example, stimulants, such as cocaine, amphetamine, and nicotine, enhance the structural complexity of the neurons, while narcotics, like morphine, decrease it (10 ).
In this review, we explore the idea that experimentally-constrained biophysical computational models that relate structural changes of neurons to different drug conditions may be able to narrow down the possible mechanisms that may underlie the neurobiology of drug addiction. Although there are a repertoire of theoretical models of drug use and addiction, they rely on theories of how the brain solves computational, information-processing, and control problems (11 , 12 ). Some of these models link these computations with specific brain areas (12 ). However, to our knowledge, they do not consider the structural and molecular changes of the neurons and, instead, focus on drug use behavior. Recently, biophysical models have informed the structural and molecular mechanisms related to learning and memory (13 -19 ). Here, we review such models and describe how they can broaden our understanding of the mechanisms behind addiction and relapse to aid the treatments for drug addicts. We focus on dendritic spines, small protrusions from the dendrites that receive excitatory input, and hence, form the postsynaptic part of the synapse. Dendritic spines are thought to be standalone biochemical computational units (20 -22 ) and it is well-documented that they undergo structural plasticity during learning and memory formation (23 ) and during drug consumption and withdrawal (5 , 7 , 24 ). Therefore, biophysical analyses of spines may provide us with the key links between these fundamental neurological processes. This work begins by giving an overview of the structure and dynamics of dendritic spines. Then, we introduce the most studied forms of synaptic changes, namely, long-term potentiation (LTP) and long-term depression (LTD), and relate them with the modifications associated with drug addiction. After this biological background, we review different types of biophysical models and discuss how they can enhance our understanding of drug addiction and other neurological diseases. We finish by describing how these biophysical models can be extended to study memory formation and storage, and identifying the opportunities for predictive model development.

Morphology
Although dendritic spines were first described by Santiago Ramón y Cajal over a century ago, it has only recently become possible, through the development of novel imaging techniques, to investigate their structure-function relationships. The dendritic spine is a bulbous protrusion from the dendritic shaft often connected by a thin neck, creating a spine head. The postsynaptic density (PSD) is located at the tip of the spine head, close to the presynaptic terminal, and it is enriched with receptors, proteins, and signaling molecules (Fig. 1A). The quantitative characteristics of the spine morphology are broadly described by: the length and width of its neck; the volume of its head; and the surface area of the PSD. These geometric features alone have been studied in terms of correlations between the morphology and function of the spine (25 -27 ). For example, the spine head volume correlates with the area of the PSD, the number of postsynaptic receptors, the readily-releasable pool of transmitters, and the neck width (25 -27 ). Moreover, dendritic spines that connect to the same axon have similar size and head volume demonstrating that they share the same activation history (28 ).
Dendritic spines also have diverse shapes. Thus, early studies classified them into three cate-gories: 1) spines with long, thin necks and small heads (thin spines), 2) spines with thicker necks and large heads (mushroom spines), and 3) short spines without well-defined necks (stubby spines) (29 ) (Fig.1B). Experiments found that spine shapes were more distinct in adults than in young animals, which suggested that stubby spines form from the dendritic shaft when the axon contacts the membrane (shaft synapse) and then turn into thin or mushroom spines (30 ). However, the current view is that a filopodia-like structure that emerges from the dendritic shaft may connect to the axon to form a synapse (31 , 32 ) evolving to a thin spine which later matures to a mushroom spine (33 ). Moreover, thin spines are more prone to respond to synaptic activity and change their morphology accordingly, suggesting that they are learning spines that can mature to mushroomlike spines, which are more stable and thus, thought to be memory spines (33 ). While mushroom spines can persist for months, thin spines only last a few days (33 ). Dendritic spines have also been divided depending on the size of their heads (34 ). On one hand, spines with small heads (filopodia and thin spines) are motile and unstable. On the other hand, spines with large heads (mushroom and stubby spines) are stable and have stronger connections. Thus, if we see spines as memory devices, large-headed spines are "write protected": they maintain preexisting connections and prevent new memory formation, while small-headed spines are "write-enabled" because they allow the acquisition of new memories (34 ).
In addition to their unique morphological features, spines have attracted considerable attention of the signaling community because it has been demonstrated they compartmentalize calcium and that spine necks can filter membrane potentials, thereby isolating biochemical and electrical synaptic input from each other (35 ). These findings have led to the hypothesis that spines implement input-specific learning rules (25 ). The shape changes in spines, in response to synaptic input, are possible by a reconfiguration of the spine cytoskeleton, mainly composed of actin filaments; actin remodeling is triggered by a cascade of chemical reactions due to calcium influx into the spine (36 -41 ). Changes in spine shape and size not only occur in development but also in adulthood (42 ), for example, during motor skill learning (43 ) and fear learning and extinction (44 ). This dynamic change to spine morphology is known as structural plasticity (45 ). Note that due to this plasticity, the shapes of spines form a continuum instead of distinct categories (26 ).

Structural Plasticity
The connection between neurons can be strengthened or weakened in a process called synaptic plasticity, and this change in the connections has a direct impact on learning and memory formation (42 , 46 -48 ). The most studied forms of these synaptic changes are long-term potentiation (LTP), which strengthens the synapse, and long-term depression (LTD), which weakens it. In LTP, glutamate released due to the high-frequency stimulation of the presynaptic neuron is captured by the α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid receptors (AMPARs) located at the PSD of the dendritic spine (Fig. 1A). Due to this glutamate influx, the spine depolarizes and releases the Magnesium (Mg 2+ ) ion of the N-methyl-D-aspartate receptors (NMDARs), which allows the calcium (Ca 2+ ) influx into the spine and triggers a cascade of reactions that increase the number of AMPARs at the PSD and AMPARs mediated current. Such increase facilitates the subsequent uptake of glutamate, and hence, the synapse strengthens. Moreover, dendritic spines show an associated increase in size, thus, linking their function and morphology (22 , 45 -47 ). LTD also depends on the activation of NMDARs, but it is induced by low-frequency stimulation of the presynaptic cell, causing glutamate release. However, the amount of glutamate is not enough to depolarize the postsynaptic cell and remove all the Mg 2+ ions blocking the NMDARs, but only a few. Hence, the Ca 2+ influx is lower than in LTP and triggers a different set of chemical reactions resulting in a decrease in the number of AMPARs at the PSD, and spine shrinkage (48 ).
How is plasticity impacted by drugs? It has been shown that cocaine use impairs the induction of LTP in spines of medium spiny neurons (MSNs) in the NAc (7 ). These spines are depressed after exposure to cocaine, which generates silent synapses that lack AMPARs but have NMDARs (4 , 5 , 7 , 24 ). The lack of AMPARs impairs LTP because these spines cannot respond to glutamate release. Because silent synapses are highly abundant during development, it has been hypothesized that the brain returns to a more juvenile state after drug exposure, which may underlie pathological drugseeking behaviors (49 ). Interestingly, AMPARs are reinserted into these silent synapses only when a stimulus related to cocaine consumption is present (7 , 24 ). This synaptic potentiation induces drug craving, and thus, may explain the high rates of relapse (7 ). Moreover, re-silencing these synapses via optogenetic removal of AMPARs inhibits relapse-like behaviors (49 ). A recent study shows that there is a time window in which previously silent synapses that are potentiated due to exposure to cues related to drug consumption can be temporally destabilized and vulnerable to disruptions before being consolidated (becoming mature) (24 ). Moreover, this window relates to the dynamics of signaling molecules that alter the spine cytoskeleton. Preventing spine maturation during this window decreases cue-induced cocaine seeking, which can be used for anti-relapse treatments (24 ) (Fig. 1C).
However, the precise neuronal and molecular substrates encoding the dynamics of drug memories have not been fully identified (4 , 24 ). Moreover, until recently, new imaging techniques allowed for the study of the mechanisms underlying the structural and molecular changes of the synapses (42 ), but they remain incompletely understood due to the complexity of the synapse and the diverse temporal and spatial scales of the experiments and underlying events.  (24 ). The red arrow signals the window when the synapse is destabilized and treatments can take place to avoid reconsolidation. Created with BioRender.com

Signaling Pathways associated with structural plasticity
Upon LTP induction, there is an increase of the number of AMPARs in the PSD due to the Ca 2+ / calmodulin-dependent protein kinase II (CaMKII) mediated phosphorylation of AMPAR auxiliary protein stargazin, which allows extrasynaptic diffusive receptors to bind to PSD95, a PSD molecule (50 -52 ). The exocytosis of extrasynaptic AMPARs into the plasma membrane depends on the rat sarcoma virus-extracellular signal-regulated kinase (RAS-ERK) pathway, RAB-GTPase proteins, soluble N-ethylmaleimide-sensitive-factor attachment protein receptors (SNARE proteins), syntaxin 4 and 13, and the motor protein myosin V (50 , 51 , 53 , 54 ). The RAS-ERK pathway is stimulated by CaMKII upon LTP induction (54 , 55 ). Thus, CaMKII has a dual function upon LTP: it immobilizes extracellular AMPARs in the PSD and replenishes the pool of extrasynaptic receptors via exocytosis. Moreover, CaMKII is one of the most important molecules for LTP. It is activated by calcium-calmodulin, which in turn is activated by Ca 2+ influx into the spine through NMDARs and it can act as a protein switch due to its autophosphorylation property. Such a state lasts longer than calcium elevation, thus CAMKII acts as a biochemical integrator of the multiple calcium pulses during LTP induction (50 ). In addition to its involvement with AMPARs, CaMKII also affects the dendritic spine cytoskeleton by stabilizing actin filaments (56 ). Actin is a globular protein (G-actin) that forms filaments (F-actin) and it is the main component of the spine cytoskeleton. F-actin is a polar structure that continuously polymerizes G-actin at its plus end and it is depolymerized at its minus end. In the spine, there is a dynamic equilibrium between G-actin and F-actin that is modulated by actin-binding proteins (ABPs) which promote F-actin depolymerization and G-actin polymerization (57 ). Moreover, the actin cytoskeleton is restructured upon LTP by an orchestrated interplay between actin and ABPs (58 , 59 ), which allows for spine enlargement. In the basal state, F-actin is bundled by CaMKII, which stabilizes the bundle. Upon NMDARs activation, CaMKII autophosphorylates and detaches from actin filaments, which allows CaMKII to interact with other molecules and remodel the actin cytoskeleton (56 , 60 ). Furthermore, CaMKII regulates actin dynamics by signaling pathways involving the Rho family of small GTPases, such as RhoA, Rac1, or Cdc42 (56 ). The activity of these GTPases is triggered by guanine-nucleotide-exchange factors (GEFs) and suppressed by GTPase activating proteins (GAPs), and they coordinately regulate ABPs activity. For example, RhoA controls profilin II and cofilin (56 , 61 , 62 ). The former binds to G-actin and facilitates the polymerization of F-actin at its plus end. Cofilin, which is also controlled by Rac1, has a dual function depending on its concentration. At low concentrations, cofilin severes F-actin and at high concentrations, it promotes F-actin nucleation (63 ). Note that the nucleation and severing events produce more actin filaments, and the combination of both events is likely to be responsible for keeping a steady distribution of the F-actin length. Thus, cofilin regulation is critical for actin dynamics. These pathways are depicted in Figure 2.
It has been shown that cocaine decreases active RhoA in the NAc MSNs, which could lead to a decrease in cofilin (64 ), affecting the spine cytoskeleton. Moreover, activation of CaMKII by D1like dopamine receptors in the NAc reinstates cocaine-seeking behavior; this is associated with an increase of AMPARs in the membrane (65 ). CaMKII also phosphorylates ∆Fos and it is required for the cocaine-mediated accumulation of ∆Fos in NAc (66 ) (red lines in Fig. 2). ∆Fos is a Fos family transcription factor that shows long-lasting accumulation in the NAc after chronic administration of any drug of abuse, thereby supporting the view that changes in gene expression contribute to drug addiction. Interestingly, it has been shown that ∆Fos upregulates the transcription of the CaMKII gene selectively in the NAc shell and that CaMKII stabilizes ∆Fos, allowing for a greater accumulation of ∆Fos that induces more CaMKII (66 ). However, how this positive feedback alters the spines is not fully understood, because of the complexity of the signaling pathways. Moreover, the signaling pathways activated by LTP induction have diverse temporal and spatial extents (46 ), and the coupling between the signaling pathways and spine structure is also complex (Fig.2). In order to disentangle the multiple spatial and temporal complexities involved in LTP induction and effects of drugs, we propose the theoretical models that incorporate biophysical details at various levels can shed light into the emergent behaviors of spines. These modeling approaches can elucidate the mechanisms underlying adaptation to addiction by allowing the incorporation of different experimental observations, eliminate possibilities that are not physicochemically feasible, and generate experimentally testable hypotheses.

Biophysical Models of Dendritic Spines
Addictive behaviors share common features with learning models (67 ), and hence, learning theories (68 ) can be adapted to investigate these behaviors. However, in this work, we focus on the structural and functional changes in dendritic spines that can inform relapse propensity. We first present a review of the models that investigate signaling pathways in the dendritic spine relevant to the mechanisms involved in drug addiction, particularly those pathways that result in structural LTP or LTD due to an increase or decrease of AMPARs and CaMKII dynamics. Then, we revise models that have incorporated idealized and realistic spine shapes linking the structural and functional characteristics of dendritic spines.

Signaling Models
Early models of synaptic plasticity, based on Ca +2 influx through the NMDARs, describe how different levels of Ca +2 (69 , 70 ) or induction protocols (71 ) lead to an increase or decrease in the synaptic weight, which describes the strength of the synaptic connection. Later models provide a more detailed description of the biochemical reactions leading to either LTP or LTD (72 ), the volume regulation of dentric spines upon LTP (73 ), or the change in AMPARs mediated current promoting synaptic modifications (74 ). These signaling models use ordinary differential equations (ODEs) to describe the biochemical reactions. They span from simple descriptions based on key assumptions (71 ) to detailed descriptions of the different chemical pathways (74 ). As technological advances allow for the study of specific signaling pathways, models have incorporated these findings to test synaptic plasticity hypotheses (75 ). The increase in computational power has also allowed the use of agent-based models (76 ), in which the molecules are agents that follow a set of rules, mimicking chemical reactions leading to LTP. Some of these models include the CaMKII activity (70 , 72 , 76 , 77 ) and GTPases (73 , 74 ).
These models test different hypotheses based on the available experimental findings and create theories about how long-term information can be stored in the brain. For example, early models test the idea that CaMKII provides a bistable switch (69 ) where "off" CaMKII is turned "on" when there is a synaptic event and stays on for long periods due to the autophosphorylation property of the kinase. Furthermore, this hypothesized switch can stay on for longer times than protein turnover because the newly synthesized unphosphorylated CaMKII subunit could become phosphorylated (78 ). The mathematical models show the feasibility of this hypothesis (69 ), and expand it to add a phosphatase switch involved in LTD, providing a tristable system to examine how a synapse can be bidirectionally modified (72 ). Models have also been used to systematically investigate how the interaction between CaMKII and other molecules can control the sensitivity of synapses to calcium signals (73 , 76 ). Moreover, the incorporation of actin remodeling revealed that it added robustness to the dendritic spine response upon LTP (73 ). One common difficulty for these models is parameter fitting since not all the reaction rates or protein concentrations are known in the spines, due to experimental difficulties. Thus, informed assumptions and decisions are made. One common assumption is that reaction rates are the same across cell types, and protein concentrations are cell-type dependent; and hence, protein concentrations are fitted in the models (74 ). Although parameters are fitted when unknown and varied to see their effects in the model dynamics, they have to be within physiological regimes (73 , 74 ).
As experiments began to investigate the source of the AMPARs that incorporate into the PSD upon LTP, i.e., whether they are exocytosed or laterally diffused into the PSD (51 , 79 -81 ), theoretical studies examined the implications of the different sources to plasticity (15 , 82 -85 ). Moreover, models started to explore how synapses are maintained over time periods that are longer than the average lifetime of these receptors and other molecules (86 -88 ), and aided the study of AMPARs clusters formation (88 , 89 ). These models propose different hypotheses for memory maintenance. Instead of relying on a bistable molecular switch, like that of CaMKII, a model of clusters of AM-PARs show that they are metastable with longer lifetimes than that of a single receptor, and that synaptic weights, which depend on the number of AMPARs, can be bidirectionally changed (88 ). An alternative hypothesis is that AMPARs exhibit a self-sustained switch related to their movement to and from the membrane (86 ). In this model, the presence of receptors at the synapse enhances the recruitment of more AMPARs, resembling the conversion from silent to active synapses (86 ). Moreover, a model based on simplified biophysical assumptions shows cooperative binding and unbinding of proteins at the PSD, replicating experimental observations (87 ). By performing model simulations, the interplay between AMPARs and CaMKII dynamics has been examined (85 , 86 ).
Another important player of structural plasticity, the actin filaments dynamics (36 , 38 ) and its dynamic equilibrium (57 ), has also been investigated providing a quantitative description of the change of the filamentous to globular actin ratio upon LTP observed experimentally (90 ). Moreover, computer simulations have aid the understanding of the F-actin binding and unbinding to CaMKII (60 ).

Models Incorporating Spine Shape
Given the structural and morphological complexity of spine shapes, the incorporation of geometrical detail is a critical aspect of computational modeling for LTP. Models started examining the role of the spatial localization of signaling molecules in synaptic plasticity (19 , 91 , 92 ), the effect of their diffusion on nearby spines (18 , 92 , 93 ) or the localization of other Ca 2+ sources, for example the endoplasmatic reticulum or the mitochondria (17 , 18 , 28 , 94 ), or astrocytes, a glia type cell of the central nervous system that contacts the synapse (95 ). Many of these models use idealized, simplified geometries where the spine head is represented by an ellipsoid or a sphere and the spine neck by a cylinder connecting the spine head with the dendrite (17 , 18 , 75 , 91 , 92 ). The dynamics of the signaling molecules are described by reaction-diffusion models that use partial differential equations (PDEs). Because some molecules are scarce in the dendritic spines, they cannot be described by their concentration amount, as is done for deterministic reaction-diffusion models, and thus, stochastic simulations were implemented to better represent their dynamics (91 , 92 ). Also, the localization and diffusion of AMPARs were examined in idealized geometries (83 , 96 ). Recent models have used reconstructed dendritic segments with spines from electron micrographs to give an insight into how the shape of dendritic spines influences calcium dynamics (28 , 97 , 98 ).
Overall, the incorporation of spatial detail in the models allowed examination of how molecular organization in the dendritic spines affect the response of signaling molecules to calcium (19 ), and the function of proteins (91 ). Moreover, different hypotheses, that are not possible to experimentally test yet, about how organelle size and contact between organelles affect Ca 2+ dynamics have been tested (17 , 18 ). Although including the spatial component in the models increases their complexity and the computational burden of the simulations, it is necessary to shed light on mechanisms that rely on the protein distribution and translocation, as structural plasticity.
Models started to explore the role of actin cytoskeleton remodeling into spine maturation, how different configurations of actin filaments derive in the different types of spines due to a balance between a force generated by actin polymerization and the resistance offered by the spine membrane (16 , 99 ), and show the importance of the spine-neck constriction by structures like actin rings (100 ) for spine stability. These models assumed isotropic forces generated by actin pushing the lipid membrane forward, which is represented as a thin elastic shell that minimizes its bending energy (101 ), to investigate the equilibrium configurations of symmetric spine shapes.
However, models that incorporate spatio-temporal dynamics of signaling often neglect the dynamics of actin remodeling because of the inherently different types of governing equations. Signaling models usually focus on reaction and diffusion in a spatio-temporal manner (17 -19 , 91 , 92 ). Models incorporating actin remodeling usually involve force balances and movement of the boundary (spine membrane in this case) in response to those forces (16 , 99 ). A model coupling the signaling spatio-temporal dynamics and the membrane dynamics resulting from a force imbalance is hard to analyze because it is not always clear how its steady state is defined: the protein location affects the spine shape that, in turn, affects the protein dynamics. Moreover, in basal conditions, proteins are not steady, they diffuse and react with other molecules, and the spine shape fluctuates (102 , 103 ). This increases the complexity of computational simulations because, at each time step, the chemical reactions and spine shape have to be updated. Nonetheless, some efforts have been made to push the boundaries on these biophysical models forward. New models were proposed to investigate these spontaneous spine shape fluctuations (104 , 105 ) and their implications in LTP (14 , 105 ) by continuously changing the membrane according to the actin dynamics. These models couple the biochemical and mechanical properties of dendritic spines, and thus, their function and structure. This allows the comparison with time-lapse microscopy data and shed light into the underlying mechanisms for spine maintenance and changes upon LTP. Moreover, it promotes testing different hypothesis regarding to the most efficient mechanisms behind synaptic function.
The level of detail also varies across studies, from simple models based on a few assumptions (71 ) to detailed models accounting for all the signaling pathways (74 ). Most of the models look at the reaction induced by Ca 2+ influx through NMDARs, but also other sources of calcium such as organelles (17 , 18 ). Moreover, the spatial representation of the dendritic spine also expands from a volume (73 ) or number of AMPARs readout (74 ) to the embodiment in static idealized (16 -18 , 28 , 83 ) and real geometries (97 ) or motile shapes (14 , 104 , 105 ).
The dynamics are described in different ways depending on the studied phenomena. To select which type of model is the most suited to investigate a biological phenomenon, it is critical to know the readout of the model and whether it is comparable with experiments. If there are no available experiments, at least there has to be evidence that such experiments can be performed with the current technology. This way, the model predictions can be tested, and the model can be used to motivate more experiments. Most signaling pathways use ODEs to describe chemical reactions (75 ) or PDEs when accounting for molecular diffusion (75 ). When the number of the studied molecules inside the spine is low, stochastic ODEs or PDEs are used instead (91 , 92 ). Other models take molecules (76 ) or filaments (90 ) as agents that follow experimentally verified rules that mimic their reactions.
Overall, there is always a trade-off between the level of detail, the spatial and temporal scale, and computational burden when proposing theoretical models. For example, for modeling the relapse propensity that lasts from hours to days (24 ), most of the fast molecular dynamics can be assumed to reach a semi-steady state. Hence, the model would only account for deviations of this state due to drug exposure. In this top-down approach, only chemical pathways affected by drug consumption are considered, and the rest are assumed to be unaffected. The affected pathways involve cytoskeleton proteins which induce remodeling. Thus, not only the variation of the AMPARs should be included, but the spine size and shape changes. Moreover, mechanical properties, such as surface area tension that regulates cell functions and triggers signals (106 ), should be considered in the models.
4 Perspective: Can we use models to enhance our understanding of synaptic changes during disease progression and to study the effect of drugs?
We consider that biophysical models can act as a computational microscope that incorporates diverse experimental findings to perform in silico experiments and test hypotheses that are not yet accessible experimentally. Therefore, they provide a platform to study the molecular effects of drug consumption in dendritic spines and simulate the impact of relapse treatments that inhibits or enhance certain signaling pathways. Although currently there is vast experimental measurements that can be use to restrict these models, these are sensitive to the the different experimental setups of drug administration: whether the drug is self-administrated or administrated by an experimenter (4 , 5 ). Thus, the models should also consider these discrepancies. Models should be suited to incorporate different spatial and temporal extents of the signal transduction underlying LTP (46 ) that is affected by cocaine. As a modeling community, this can be achieved by building modular, shareable models for different time-and spatial scales that can be coupled. By developing biophysical models to understand the changes in the synapse during drug abuse, we can gain insight into mechanisms that may be involved in other neurological diseases. For example, dendritic spine alterations have been found in autism spectrum disorders, schizophrenia, and Alzheimer's disease (49 , 107 ). Such alterations correlate with the expression of proteins associated with the cytoskeleton (108 ), or the number of AMPARs (109 ). The models presented here can be extended to study the complex interplay between membrane mechanics, cytoskeleton dynamics, and protein interaction during endo-and exocytosis of extracellular vesicles (EVs) containing AMPARs that are inserted and removed from the spine membrane, respectively. Although such mechanisms have been modeled for other types of cells (110 -112 ), they have not been incorporated into dendritic spine studies. Furthermore, EVs containing proteins and RNAs are also endo-and exocytosed. Recent evidence shows EVs are involved in pathologies of neurological diseases (113 , 114 ), and that they can potentially serve as diagnostic tools (113 -115 ) or as vehicles for medicine delivery (116 ). Therefore, biophysical models can be used to complement the efforts in understanding the underlying mechanisms to propose new avenues to diagnose and treat neurological diseases.

Outlook
We have assumed that dendritic spines are capable of storing memories: they protrude from dendrites and then connect to presynaptic sites. After LTP, the number of AMPARs at the PSD increase, enhancing the AMPAR-mediated conductivity and enlarging the spine. These mature spines can last over days or disappear to promote the creation of new memories, allowing for dynamic memory (117 ). But the idea that memories reside on synapses has not been conclusively proven yet (118 ). Moreover, synaptic proteins have half-lives of a few days and there is a continuous turnover of spines (118 , 119 ), but memories can last years. This leads to an alternative hypothesis stating that memories are rather stored in engrams, that are populations of neurons activated by an experience, across different brain regions (120 ). There have been modeling efforts to understand these different mechanisms, from the study of different learning rules that dictate the strength changes of the connection between neurons (68 , 121 ) to investigate how memories can be stored and recalled in different cell assemblies (122 ) or how these networks can learn without disrupting memories (123 ). Besides spanning across the spatial scale of the synapse and cell assemblies, the mechanisms underlying learning and memory storage span across different time scales (124 ). Although models are addressing how memories can be stored despite high spine turnover (13 ) or how multiple spines in the dendrite contribute to single neuron computations and how it is affected by dendritic morphology (125 , 126 ), more studies are needed to address the pressing problems facing society.