Transcriptional correlates of memory maintenance following long-term sensitization of Aplysia californica.

We characterized the transcriptional response accompanying maintenance of long-term sensitization (LTS) memory in the pleural ganglia of Aplysia californica using microarray (N = 8) and qPCR (N = 11 additional samples). We found that 24 h after memory induction there is strong regulation of 1198 transcripts (748 up and 450 down) in a pattern that is almost completely distinct from what is observed during memory encoding (1 h after training). There is widespread up-regulation of transcripts related to all levels of protein production, from transcription (e.g., subunits of transcription initiation factors) to translation (e.g., subunits of eIF1, eIF2, eIF3, eIF4, eIF5, and eIF2B) to activation of components of the unfolded protein response (e.g., CREB3/Luman, BiP, AATF). In addition, there are widespread changes in transcripts related to cytoskeleton function, synaptic targeting, synaptic function, neurotransmitter regulation, and neuronal signaling. Many of the transcripts identified have previously been linked to memory and plasticity (e.g., Egr, menin, TOB1, IGF2 mRNA binding protein 1/ZBP-1), though the majority are novel and/or uncharacterized. Interestingly, there is regulation that could contribute to metaplasticity potentially opposing or even eroding LTS memory (down-regulation of adenylate cyclase and a putative serotonin receptor, up-regulation of FMRFa and a FMRFa receptor). This study reveals that maintenance of a "simple" nonassociative memory is accompanied by an astonishingly complex transcriptional response.

Long-term memories persist despite continuous molecular turnover. This seems possible, in part, due to sustained structural plasticity (Caroni et al. 2012). How is this growth maintained? A number of mechanisms have been proposed, but one requirement seems to be transcriptional regulation that extends beyond memory induction. Specifically, the encoding of long-term memories is accompanied by multiple waves of changes in gene expression (Barzilai et al. 1989; Alberini 2009). Moreover, work in a number of model systems has shown that the maintenance of long-term memory can be impaired by blocking transcription during critical periods after training (Igaz et al. 2002; Lefer et al. 2012). What are the targets for maintenance-related changes in gene expression? Surprisingly, this has not yet been fully elucidated in some memory paradigms (see below) or even for common forms of long-lasting synaptic plasticity (Abraham and Williams 2008).

Here we characterize the transcriptional correlates of memory maintenance following long-term sensitization (LTS) training in the marine mollusk Aplysia californica. Aplysia have long served as an attractive model organism for studying the molecular mechanisms of memory. One particular focus has been LTS (Pinsker et al. 1973), a learning paradigm in which repeated exposure to a noxious stimulus produces a long-lasting, transcription-dependent increase in reflex responsiveness (Castellucci et al. 1989). Sensitization of the tail-elicited siphon-withdrawal reflex provides an especially attractive system for transcriptional analysis because (1) sensitization can be applied and expressed unilaterally (Scholz and Byrne 1987), allowing for powerful within-subjects comparisons, (2) sensitization memory is known to depend at least in part on physiological changes in the VC nociceptors of the pleural ganglia (Cleary et al. 1998), providing a behaviorally relevant target for transcriptional analysis, and (3) transcriptional and behavioral changes can be correlated at the level of individual animals, allowing exploration of individual differences in memory retention (Bonnick et al. 2012).

We have previously analyzed the transcriptional correlates of LTS encoding (Herdegen et al. 2014b), showing that 1 h after training there is a strong up-regulation of 81 transcripts, including those encoding transcription factors and transcription-factor regulators: ApC/EBP (GenBank: U00994; Alberini et al. 1994), ApCREB1 (GenBank: NM_001256437; Bartsch et al. 1998), ApEgr (GenBank: KC608221; Cyriac et al. 2013), and ApC/EBPγ (GenBank: EB233406).

How might the transcriptional response accompanying encoding change during the maintenance of an LTS memory? This remains unclear, but there have been several previous proteomic screens. In one, LTS was induced over a 4 d period and abdominal ganglia were harvested 1 d after training (Castellucci et al. 1988); this led to the identification of four reliably regulated proteins. These same two proteins were reconfirmed as persistently regulated in VC sensory neurons in a second screen that used serotonin exposure to mimic the induction of LTS memory (Barzilai et al. 1989). Finally, a third screen was recently conducted examining changes in protein expression 1 and 2 d after serotonin exposure in cultured Aplysia ganglia; this identified 18 up-regulated and 1 down-regulated proteins 1 d after exposure in the pleural ganglia (Monje et al. 2012).

These screens plus additional gene-of-interest studies have identified several changes at the level of mRNA expression that accompany the maintenance of LTS memory. Some changes persist from the encoding phase (Herdegen et al. 2014a), including up-regulation of ApEgr (GenBank: KC608221; Cyriac et al. 2013), ApGlyT2 (sodium- and chloride-dependent glycine-dependent transport 2, GenBank: XM_005092349), ApVPS36 (Vacuolar protein-sorting-associated protein 36-like) and an uncharacterized transcript (LOC101862095, GenBank: XM_005113453). Others emerge after encoding, including a delayed increase in the expression of ApBiP (GenBank: NM_001204652; Kuhl et al. 1992), ApCalreticulin (GenBank: NM_001204594; Kennedy et al. 1992), ApTBL-1 (GenBank: NM_001204563; Liu et al. 1997), and Sensorin (GenBank: NM_001204654; Schacher et al. 2000). There are probably additional transcriptional changes. LTS training can induce tremendous outgrowth of Aplysia sensory neurons (Bailey and Chen 1983; Wainwright et al. 2004); this likely requires the regulation of many transcripts. In addition, the maintenance of LTS memory requires persistent changes in DNA methylation, as temporary inhibition of DNA methyltransferase can eliminate the expression of LTS memory even days after induction (Pearce et al. 2017). This indicates the requirement for methylation-mediated changes in gene expression for the maintenance of LTS, though the targets remain unclear.

To fully characterize the transcriptional correlates of LTS sensitization we used microarray analysis to measure the changes in gene expression persisting 24 h after training. We analyzed changes in the pleural ganglia which contain the VC nociceptors (Walters et al. 1983) and which are thought to mediate much of the expression of LTS memory (e.g., Walters 1987a) as well as interneurons contributing to withdrawal circuitry. Analyzing the whole ganglia reflects contributions related to both generalized- and site-specific sensitization (see Materials and Methods). As microarray analysis is exploratory, we took several steps to ensure results would be reproducible: preregistration of our microarray analysis plan, sample size planning for adequate power, and qPCR validation in an independent sample. We find that the transcriptional correlates of LTS are remarkably complex, involving regulation of an appreciable fraction of the Aplysia genome.

Results

LTS training produces unilateral LTS that lasts >4 d

To induce a LTS memory, animals received a series of noxious shocks to one side of the body (Fig. 1). As previously reported, this led to a unilateral LTS memory (Fig. 2). On the trained side, T-SWR durations averaged 8.9 sec prior to training and 16.5 sec when measured 24 h after training, an average increase of 7.5 sec (95% CI [6.9, 8.2], dunbiased = 4.1 95% CI[3.3, 5.0], r = 0.32, t(49) = 24.4, P < 0.001). On the untrained sides, SWR reflexes were 9.0 sec prior to training and 8.7 sec when measured 24 h after training, a slight decrease of 0.3 sec (95% CI[−0.7, 0.2], dunbiased = −0.2 95% CI[−0.6, 0.2], r = 0.12, t(49) = −1.03, P = 0.31). Thus, we observed the expected interaction between side of training and phase of testing (MTrainedDiff-UntrainedDiff = 7.8 sec 95% CI[7.1, 8.5], dunbiased = 3.9, 95%CI [3.1, 4.8], r = 0.20, t(49) = 22.1, P < 0.001). All 50 animals met our preregistered quality control of showing at least a 30% increase in T-SWR duration on the trained side.

Figure 1.

Long-term sensitization of the tail-elicited siphon-withdrawal reflex (T-SWR). (A) Cartoon diagram of the body of an Aplysia. T-SWRs are evoked by applying an innocuous shock to the left or right tail (arrows). The duration of the T-SWR serves as an index of behavioral responsiveness. For LTS training, a noxious shock is applied along the length of one side of the body (lightning bolts). (B) Experimental protocol. First, baseline T-SWR measures are made on the left and right side of the tail, then LTS training is applied to one side of the body, then T-SWR measures are made again 24 h after training. Immediately after post-tests, pleural ganglia from the trained and untrained side are harvested in matched pairs of left-trained and right-trained animals. The pleural ganglia contain the somas of the VC nociceptors which help mediate sensitization memory.

Figure 2.

LTS training produces unilateral sensitization memory. T-SWR durations before (Baseline) and 24 h after (24 h Post-Test) LTS training on both the trained (left panel) and untrained sides (right panel). Individual data from each animal is shown with open circles with a line joining that animal's Baseline and 24 h Post-Test responses. Averages from each testing phase are shown as solid circles with 95% confidence intervals. Each difference score from Baseline to 24 h Post-Test is shown with an open triangle; solid triangles represent average differences scores with 95% confidence intervals.

To confirm that our analysis represents the maintenance phase of memory, we also completed a time-course analysis in an independent set of animals. Animals received the same 1-d LTS training but were tested for retention 1, 2, 4, 7, and 9 d after training (n = 15). As shown in Figure 3, our 1 d LTS protocol produces a unilateral memory that is strongly expressed at 1, 2, and 4 d after training (dunbiased = 3.23, 2.6, and 3.3, respectively; each interaction test significant at P < 0.001). At 7-d, the evidence for remaining memory reached statistical significance (MTrainedDiff-UntrainedDiff = 2.1 sec 95% CI[0.35, 3.9], dunbiased = 0.8, 95%CI [0.1, 1.5], r = 0.21, t(14) = 2.6, P = 0.02). Notably, the effect size was much smaller at 7 d, and there were large individual differences in retention. At 9 d, there was no longer strong evidence for LTS memory (MTrainedDiff-UntrainedDiff = 0.7 sec 95% CI[−0.5, 1.9], dunbiased = 0.3, 95%CI [−0.2, 0.9], r = 0.41, t(14) = 1.2, P = 0.24).

Figure 3.

LTS Memory persists for >4 d after training. T-SWR durations are shown normalized as a percentage of baseline 1, 2, 4, 7, and 9 d after LTS training. Black squares and open circles with error bars represent mean scores and 95% CI for trained and untrained sides, respectively. The semitransparent lines with small squares and circles represent individual data points.

LTS training increases expression of ApBiP and ApEgr

LTS training produces a delayed but long-lasting increases in the expression of ApBiP (Kuhl et al. 1992) as well as a rapid and long-lasting increase in the expression of ApEgr (Cyriac et al. 2013). As a quality control, we used qPCR to confirm that training had produced these expected transcriptional responses. Tissue was harvested immediately after the 24 h post-tests. Samples from two animals trained on opposite sides were pooled (50 animals → 25 sets of samples). Issues with tissue processing caused 3 sets to be discarded, leaving 22 for analysis.

As expected, expression of ApBiP was strongly regulated by LTS training (Fig. 4), with a mean fold change (MFC) of 1.9, indicating nearly double the expression on the trained side relative to the untrained side (95% CI[1.5, 2.3], dunbiased = 1.0 95% CI[1.2, 2.6], t(21) = 6.18, P < 0.001). Similarly, ApEgr was also strongly regulated, with nearly three times the expression on the trained side relative to the untrained side (MFC = 2.9 95% CI[2.1, 4.1], dunbiased = 2.0 95% CI[1.2, 2.7], t(21) = 6.54, P < 0.001).

Figure 4.

LTS training up-regulated the expression of ApBiP and ApEgr. Fold-change scores (trained/untrained) in the expression of ApBip (left) and ApEgr. Open symbols represent individual samples; solid symbols are group means with 95% confidence intervals. The dashed line represents a fold-change score of 1, meaning equal expression on the trained and untrained sides. Note the use of a log-scale axis to give equal visual weight to up- and down-regulation.

Although training produced strong regulation of ApBiP and ApEgr, three sets did not have higher ApBiP and ApEgr expression on the trained side. Following our preregistered quality controls (Herdegen et al. 2014a), these samples were not further analyzed, leaving 19 sets. Of these, 8 were used for microarray analysis; the other 11 were held back for independent confirmation with qPCR.

The transcriptional correlates of LTS maintenance are complex

To elucidate the transcriptional correlates of memory maintenance, we conducted microarray analysis on eight samples. Each array was conducted with a two-color approach contrasting paired trained and untrained samples. To ensure practical as well as statistical significance, transcripts were marked as regulated only when there was clear evidence of more than a 10% change in expression in either direction. Correction for multiple comparisons was also used to limit the false-discovery rate (FDR) to 5%.

From the microarray analysis 1198 transcripts were identified as strongly regulated (4.6% of the 26,091 unique transcripts tested; Supplemental Table 1). Of these, 748 were up-regulated; 450 were down-regulated.

Following best practices (Allison et al. 2006), we used qPCR in an independent sample to validate the microarray results. Specifically, we measured the expression of 43 different transcripts in the remaining 11 sets of samples. Overall, we found very strong convergent validity (Fig. 5, r = 0.89 95% CI[.81, 0.94], N = 43, P < 0.001). Of 30 transcripts flagged as significantly regulated in the microarray all 30 were significantly regulated using qPCR in the independent sample (false positive rate = 0/30 = 0%, Supplemental Fig. 1).

Figure 5.

Predictive validity of the microarray analysis. This scatter plot shows the fold-change scores (trained to untrained) for 44 transcripts measured via microarray (y-axis, N = 8) and in an independent sample using quantitative PCR (x-axis, N = 11). The intersection of the two axes is 1, which represents no change in gene expression.

Might the list of regulated transcripts be incomplete? To address this issue we estimated the likely proportion of true negatives in the microarray data. We used the approach developed by Langaas et al. (2005) which is based on analysis of the distribution of a large set of p values. This analysis indicated that 94.5% of tested transcripts are likely to be truly nonregulated. Given this, the false-negative rate is estimated to be 0.9% (100%—4.6% regulated—94.5% true negative); this is equivalent to an additional 237 transcripts. Thus, our list is probably not complete, but it is likely our analysis identifies the majority of strongly regulated transcripts represented on the array.

Up- and down-regulated transcripts are associated with different biological functions

To help make sense of the list of regulated transcripts we used gene ontology analysis. Using standard settings in Blast2Go yielded a relatively meager set of GO annotations: only 274 of 748 (37%) up-regulated transcripts and 80 of 450 (18%) down-regulated transcripts. Therefore, the analyses reported here must be considered provisional.

For up-regulated transcripts, the most common biological processes regulated were transport, translation, signal transduction, and a number of terms related to the regulation of gene expression (Fig. 6A). For down-regulated transcripts, the most common biological processes GO terms were catabolic process, cellular amino acid metabolic production, intracellular signal transduction, and cofactor metabolic process (Fig. 6B). Examples of transcripts annotated in these categories are provided in Supplemental Table 2.

Figure 6.

Most-represented biological process gene ontology for regulated transcripts. Most represented terms are shown for up-regulated transcripts (Panel A) and down-regulated transcripts (B). These graphs show the most-represented terms across gene ontology levels (multilevel); the numbers in parenthesis represent a scoring/weighting system used by Blast2Go (Conesa et al. 2005) to weight the representation of gene ontology terms.

Comparing the biological process terms in these two sets (overall FDR ≤ 5%) showed that up-regulated transcripts were significantly enriched for a number of terms related to regulating the production of proteins (Table 1). This included basic regulation of gene expression, regulation of RNA processing, translation, and several other GO terms related to macromolecule/peptide synthesis and processing. Down-regulated transcripts showed an enrichment for transcripts related to catabolism and the regulation of metabolism. Thus, LTS memory is associated with a widespread up-regulation of transcripts required to produce proteins and a down-regulation of protein breakdown—both factors that might be expected to promote growth.

Table 1.

Biological process terms over-represented when comparing up- and down-regulated transcripts

To supplement this analysis, we also matched microarray probes to genomic mRNAs from the current draft of the Aplysia genome. Of the 1198 regulated transcripts, 671 (52%) were definitively matched to an mRNA, and 552 of these were (46% of total) were matched to named mRNAs. We manually scanned and collated this set (Lakhina et al. 2015). Tables 2 and 3 present some of the major themes that seemed evident from this approach; the Discussion section comments on each category and provides connections to prior research.

Table 2.

Curated list of late-regulated transcripts related to the production of protein

Table 3.

Curated list of late-regulated transcripts potentially related to the expression of LTS memory

The transcriptional correlates of encoding and maintenance are almost completely distinct

To what extent are memory encoding and memory maintenance transcriptionally distinct? To answer this question, we compared the microarray results obtained 24 h after training with a previous analysis of the rapid (1 h) transcriptional response to LTS training (Herdegen et al. 2014b). In terms of regulation due to training, there was essentially no correlation in raw log-fold-change scores across time points (r = 0.03). This can be an underestimate of the true correlation, however, due to related measurement error in both data sets. When adjusting for this possibility using the genas function (“genuine association,” Ritchie et al. 2015a) in limma we found that there is a weak positive relationship in regulation from encoding to maintenance (Fig. 7A, r = 0.25, N = 26,091). This indicates that at least some of the pattern of regulation at 1 h is preserved at 24 h, but only enough to account for 6% of the variance in the 24 h data. This weak correlation was not due to poor measurement, as overall expression levels were very consistent across the two data sets (Fig. 7B, r = 0.92, N = 26,091).

Figure 7.

Relationships in Gene Expression 1 and 24 h after LTS Training. (A) Comparison of adjusted fold-change scores (trained to untrained) between time points. This panel represents the degree to which learning-induced changes in expression are correlated from encoding to maintenance. Transcripts regulated at both time points appear in the upper-right or lower-left quadrants. This panel was produced using the genas function from limma (Ritchie et al. 2015a), which adjusts fold-change scores to correct for potential measurement error. (B) Comparison of overall expression (transformed into Z-scores) between time-points. Overall expression scores represent an average across trained and untrained sides and simply indicate the overall degree of expression in the pleural ganglia. The line of best fit is for a simple linear regression.

The weak positive association across time points indicates that a few transcripts are regulated both rapidly and persistently. To identify these transcripts, we compared the lists of strongly regulated transcripts at the 1 and 24 h time points. Consistent with the low overall correlation in regulation, we found very little overlap. Nevertheless, we found that 15 of the 81 (19%) transcripts that are strongly regulated at 1 h after training continue to be strongly regulated 24 h after training.

To ensure that identification of these “rapid-and-persistent” transcripts was not incomplete due to insufficient power, we reanalyzed the 24-h microarray data using only transcripts regulated at 1 h, and vice versa. These more-focused analyses increase power by decreasing the number of comparisons being made, and thus reducing the degree of correction needed for multiple comparisons. These analyses identified an additional 11 transcripts that could be regulated both rapidly and persistently (for a total of 26, Table 4). One of these, however, was regulated in the opposite direction across time points (GenBank: EB316959.1; up-regulated at 1 h, down-regulated at 24 h).

Table 4.

Putative maintenance transcripts: list of transcripts regulated rapidly and persistently by LTS training

Notably absent from the list of transcripts regulated at both time-points was ApCREB1, which has been proposed to play a role in memory maintenance through a transcriptional positive feedback loop (Liu et al. 2008, 2011a) and which shows increased mRNA expression 24 h after serotonin treatment in isolated ganglia (Liu et al. 2008). Although the expression of ApCREB1 is sharply up-regulated 1 h after LTS training (Herdegen et al. 2014a,b), we did not observe persistent (24 h) up-regulation of ApCREB1 in the microarray data (MFC = 1.1, pcorrected = 0.99). To confirm this we also measured ApCREB1 in our validation set using qPCR. This approach also failed to show up-regulation of ApCREB1 mRNA (Fig. 8, MFC = 1.0 95% CI [0.7, 1.4], dunbiased = 0.2 [−0.7, 1.2], t(10) = −0.08, P = 0.94), which is consistent with previous findings from our laboratory (Herdegen et al. 2014a). ApCREB1 could remain elevated at the protein level (Liu et al. 2011a).

Figure 8.

ApCREB1 expression 24 h after LTS training. Data are with qPCR from the validation set. Open circles are individual samples, the filled circle is the mean with 95% CI.

We also did not observe persistent (24 h) regulation of ApCREB2 mRNA (GenBank: NM_001204701; Bartsch et al. 1995) with the microarray (MFC = 1.2, pcorrected = 0.11). In culture, serotonin treatment produces a down-regulation of ApCREB2 mRNA at 24 h in both Aplysia sensory neurons (Rajasethupathy et al. 2012) and motor neurons (Hu et al. 2015). However, it has proven difficult to detect these changes with serotonin treatment in whole ganglia (Liu et al. 2011b). We also have not observed regulation in the pleural ganglia 24 h after LTS training (MFC = 0.97 95% CI[0.6, 1.2], dunbiased = −0.24 [−0.8, 0.3], t(12) = −0.89, P = 0.39; unpublished data from Bonnick et al. 2012), though note that the CI obtained is long. It could be that different patterns of ApCREB2 regulation in different neuronal subtypes could be obscured at the level of whole ganglia.

We conclude that there is a small core of transcripts that are rapidly, persistently, and consistently regulated by LTS training. Otherwise, the transcriptional correlates of encoding and maintenance seem to be entirely distinct.

Discussion

We previously reported that the encoding phase of LTS memory (1 h after training) is associated with strong up-regulation of 81 transcripts (Herdegen et al. 2014b). We now show that the transition to memory maintenance (24 h after training) is associated with a tremendous ramification of this response, with clear evidence of strong regulation of 1198 unique transcripts. These two phases of transcriptional regulation are almost completely distinct, with very little of the encoding phase pattern of regulation preserved during maintenance.

The transcriptional regulation observed during maintenance is not only distinct; it is also remarkably complex. The current draft of the Aplysia genome encompasses 21,426 gene models. Thus, maintenance of a “simple” nonassociative memory seems to involve regulation of an appreciable fraction (5.6%) of the entire genome. This is probably not an overestimation, as we used strict controls against Type I error and a transcriptome database to de-duplicate our EST array. Rather, we have likely underestimated the degree of regulation, due to false negatives (estimated to be ∼200 additional transcripts). In addition, our array is thought to represent only 50%–60% of neuronally expressed genes, and could also fail to detect differential regulation across neuronal subtypes within the pleural ganglia. Assuming array probes are representative, LTS training could mobilize changes in the expression of >10% of all Aplysia genes

The complexity of response we observed here was not captured in previous proteomic screens, which identified only a handful of changes in protein expression 1 d after LTS training (Castellucci et al. 1988) or serotonin exposure (Barzilai et al. 1989; Monje et al. 2012). It is possible that the transcriptional changes we observed do not always yield alterations in protein expression. More likely, though, is that this reflects a difference in sensitivity across techniques. Despite this, there is a strong concordance between our results and these previous screens. For example, our list overlaps with at least 5 of the 19 proteins identified by Monje et al. (2012): calponin, calreticulin, glutathione S-transferase, a heterogeneous nuclear riboprotein, and tubulin.[1]

Given the complexity of the transcriptional correlates of memory maintenance, a central issue becomes how to make sense of it. Here we use time-course analysis to roughly partition transcripts into different functional groups.

Three hypothesized types of memory mechanisms

Most research into the mechanisms of memory is predicated on the assumption that learning leads to changes in protein synthesis, and that these changes serve as a trigger to launch the cellular and molecular programs that will then consolidate the long-term memory. According to this “consensus” model (Klann and Sweatt 2008), there are three general types of mechanisms at play: (1) triggering mechanisms initiated at encoding which spark the consolidation process, (2) core maintenance mechanisms, which are activated by the triggering mechanisms and then self-perpetuate to maintain the memory, and (3) effector mechanisms, which are sustained by core maintenance mechanisms, and which are the cellular and synaptic changes enabling the physiological and behavioral readout of the memory.

This framework suggests that the transcripts regulated after LTS training primarily serve one of these three functions (triggering, maintaining, or effecting). Although definitive assignment to one of these categories requires experimental manipulation, a heuristic categorization can be achieved through time-course analysis. Specifically, triggering mechanisms are thought to be rapid but transient, occurring at the point of encoding but lasting only long enough to initiate the startup of the core maintenance mechanisms. Effector mechanisms, on the other hand, should exhibit a delayed but persistent time course, as they are not initiated until the core maintenance mechanisms begin to be established, but persist as long as the maintenance mechanisms self-perpetuate. Finally, the core maintenance mechanisms are expected to be semirapid (following just after triggering mechanisms) and persistent. These assumptions are lucidly explicated in “two loop” models of memory maintenance (Smolen et al. 2009; Zhang et al. 2010). These models feature a fast auto-activation loop (core maintenance mechanisms) that is rapidly and persistently activated at the time of learning and a slower synthesis loop (effector mechanisms) that is slowly but persistently activated.

Using this framework, we assign “provisional” functions to transcripts based on if they are regulated only at 1 h after LTS training (putative triggering transcripts), only at 24 h after LTS training (putative effector transcripts), or at both time points (putative maintenance transcripts). Below we consider in more detail both the putative maintenance and effector transcripts.

Rapid and persistently regulated transcripts: core maintenance mechanisms?

We identified only 25 transcripts that are regulated in a manner suggestive of maintenance functions (rapid, persistent, and consistent; Table 4). Supporting this classification, two transcripts on this list have already been strongly implicated in memory maintenance.

One of these is ApEgr (aka Zif268/zenk/TIS8/krox-24/NGFI-A). Egr is a transcription factor that is essential for long-term potentiation (LTP), several forms of memory in mammals (for a review, see Poirier et al. 2008), and song-learning in avian species (Moorman et al. 2011). The degree to which Egr is up-regulated after learning is correlated with the expression of both LTP in mammals (Abraham et al. 1993) and LTS memory in Aplysia (Cyriac et al. 2013).

A second potential maintenance gene is a putative Aplysia homolog of TOB1 (transducer of ErbB-2; GenBank: EB349048). TOB1 is necessary for some forms of memory in rats (Jin et al. 2005; Wang et al. 2006), and is a key regulator of CPEB (Hosoda et al. 2011; Ogami et al. 2014). In Aplysia, CPEB functions as a core maintenance mechanism at individual synapses, as it is switched into a self-perpetuating active state during long-term facilitation and plays a critical role in supporting the local protein synthesis required for the persistence of long-term facilitation (Si et al. 2003; Miniaci et al. 2008). We did not observe transcriptional regulation of CPEB1 nor of a CPEB2-like transcript (our array probe for CPEB1 does not distinguish between its isoforms).

Most of the remaining rapid-and-persistent transcripts are uncharacterized, providing exciting avenues for future research. However, functional assignment by the time course of regulation is heuristic and likely to be wrong in some cases. For example, the rapid-and-persistent up-regulation of a FMRF-amide transcript (FMRFa, GenBank: M11283) seems unlikely to be related to memory maintenance; this is discussed further below. Another important caveat is that the core maintenance mechanisms may not be purely transcriptional; there may also be persistent mechanisms related to both post-translational (Si et al. 2003) and epigenetic (Pearce et al. 2017) mechanisms.

Delayed but persistent transcripts: what are the effectors affecting?

The vast majority of the transcripts we identified in this screen 24 h after LTS training were not strongly regulated in our previous screen 1 h after LTS training (Herdegen et al. 2014b). We thus provisionally characterize these transcripts as effectors. So what exactly, is being affected?

Increase in the machinery of protein production

Although making sense of such a large set of transcripts is daunting, there is a clear and widespread up-regulation of the machinery of protein production (see Table 2). This spans all levels of the central dogma, with up-regulation of transcripts related to transcription, RNA processing, translation, folding, packaging, and transport. This makes sense given the enormous outgrowth of sensory neurons that can be induced by LTS training (Bailey and Chen 1983, 1989; Wainwright et al. 2002). Notably, this up-regulation includes transcripts encoding portions of several eukaryotic translation initiation factors (eIF1, eIF2, eIF3, eIF4, and eIF5) as well as subunits of eIF2B, the guanine nucleotide exchange factor (GEF) critical for eIF2 function. It is already well-appreciated that regulation of translation is critical for long-term memory (for review, see Klann et al. 2004; Costa-Mattioli et al. 2009; Wang et al. 2010). For example, eIF3e facilitates the induction of LTS in Aplysia kurodai (Lee et al. 2008),[2] eIF4 is critical for nociceptive plasticity in mice (Melemedjian et al. 2010), and enhancing eIF2B function enhances memory formation in rats (Sekine et al. 2015). What is notable about our data is how widespread and consistent this up-regulation is; it does not seem highly targeted.

The up-regulation in transcripts related to protein production did not seem to extend to the ribosomes. The array contained 58 transcripts annotated as putative ribosomal proteins and 2 annotated as rRNAs; none of these were flagged as strongly regulated. There was up-regulation, though of transcripts annotated as ribosomal binding protein 1-like (GenBank: EB274483) and ribosomal protein S6 kinase δ-1-like (GenBank: EB252768). It may be that LTS training does not induce regulation of ribosomal transcripts, or that the regulation is at a different time point, or that it is too subtle to easily detect.

One important caveat is that the regulation of translation can be distinctive in different cell types in the pleural ganglion (Dyer and Sossin 2013), so this list of transcripts could represent the union of several distinctive pathways regulated by LTS training in the samples we collected. We have just begun to explore this issue by conducting additional microarray analyses on isolated VC sensory neuron somata. Although still at a pilot stage (N = 4), we have so far found that regulation in the VCs is very similar to what is observed in the pleural ganglia as a whole (r = 0.79 estimated genuine association, N = 26,091, data not shown).

Activation of the components of the unfolded protein response

Consistent with a global increase in protein production, there also seems to be transcriptional up-regulation of several components of the unfolded protein response (UPR). This includes the previously discovered strong up-regulation of BiP/GRP78 (Kuhl et al. 1992) as well as the new discovery of up-regulation of transcripts encoding putative homologs of CREB3/Luman (CREB-3 like protein, GenBank: EB252287), CREB3 regulatory factor (GenBank: EB297571), eIF-2-α kinase activator GCN1 (GenBank: EB230807), and AATF3 (apoptosis-antagonizing transcription factor 3, GenBank: EB259679). The UPR maintains homeostasis of protein production, and initial activation produces an increase in endoplasmic reticulum function (for review, see Walter and Ron 2011). In the CNS, there is growing appreciation for the important role the UPR plays in diverse neuronal functions (Godin et al. 2016). For example, in rats peripheral nerve injury activates the UPR in the cell bodies and axons of injured sensory neurons; blocking this response impairs regenerative axonal outgrowth (Ying et al. 2015). It is difficult to assess how widespread this potential UPR activation is, as many key transcripts are not yet fully characterized in Aplysia.

Changes in neuronal function

What are all the new proteins being produced? There seem to be changes across diverse sets of transcripts, many of which seem likely to help mediate the expression of LTS memory (see Table 3). For example, we observed up-regulation of many cytoskeleton-related transcripts, which could be important for the structural plasticity that can occur following LTS training (Bailey and Chen 1983). We also observed mixed regulation of several transcripts related to pathfinding and synaptic targeting: down-regulation of a putative homologs of fasciclin (GenBank: EB225102) and headcase (GenBank: EB255051), and up-regulation of two semaphorin-like transcripts (GenBank: EB252581, EB276728). In addition, there was strong regulation of transcripts likely to encode key synaptic proteins, transmitter and growth factor receptors, transmitter transporters, and ion channels.

Although the transcriptional events we observed are global, they could have synapse-specific effects if the mRNAs are targeted for local translation at synapses. Consistent with this possibility, we observed regulation of a number of transport-related transcripts. This includes expression changes in several transcripts that seem to encode heterogeneous nuclear ribonucleoproteins (GenBank: EB250478.1, EB314164.1, GD227870.1, EB228893.1) and one that encodes a KIF protein (GenBank: EB225867.1). As these proteins are all implicated in neuronal mRNA transport (Hirokawa 2006), regulation of their expression could relate to targeting translational changes to specific synapses (Wang et al. 2010).

Supporting the validity of our approach, regulated transcripts include nearly all those previously identified as persistently regulated following LTS training (marked with an * in Tables 2–4) as well as a number of transcripts previously linked to long-term memory or long-term neuronal plasticity in other systems. For example, we observed up-regulation of transcripts encoding putative homologs of menin/SCG2 (menin like, GenBank: EB264888.1) IGF2 (insulin-like growth factor 2, GenBank: GD220814), and FGFr-4 (fibroblast growth factor receptor 4-like, GenBank: CK323234). Menin up-regulation is essential for spinal plasticity expressed after peripheral nerve injury (Xu et al. 2012). Insulin-like growth factor 2 (IGF2) helps mediate inhibitory-avoidance learning in rat and can even serve as a memory enhancer (Chen et al. 2011). Some fibroblast growth factor receptors are essential for memory maintenance and LTP (Zhao et al. 2007). The FGFr we observed regulated is known not to influence serotonin-induced translocation of PKC (GenBank: XM_013087075; Farah et al. 2016) but may play other roles in LTS memory

Finally, our list contains a large number of novel transcripts, many of which make sense relative to what is already known about the mechanisms mediating LTS memory. For example, we observed strong down-regulation of a transcript that seems to encode a specific voltage-gated potassium channel subunit. This could help mediate the persistent spike narrowing which accompanies LTS training (Antzoulatos and Byrne 2007).

Metaplasticity

Another noticeable pattern in our results is the down-regulation of some transcripts known to be associated with the induction of LTS memory. Specifically, there seems to be a substantial down-regulation of three of the four known Aplysia homologs of adenylate cyclase (AC-AplA, GenBank: AY843027; AC-AplC, GenBank: AY843026; AC-AplD, GenBank: HM030824). The A isoform is expressed in sensory neurons and is activated by bath application of serotonin (Lin et al. 2010), a treatment which mimics the induction of LTS memory. There is also a down-regulation in a transcript encoding a type-4 cAMP-specific phosphodiesterase (PDE4, Genbank: FF071603.1). PDE4 limits the function of cAMP, so down-regulation might be expected to foster the induction of LTS. However, work in Aplysia kurodai found that both over- and underexpression of PDE4 limits the ability of serotonin to activate PKA and blocks serotonin-induced long-term facilitation (Park et al. 2005). Finally, we also observed a decrease in the expression of a putative serotonin-receptor transcript (GenBank: EB255570.1), although it is unclear if this is one of the receptors required for induction of LTS. Overall, it seems as through these changes represent a form of metaplasticity—a change in the ability to induce or express plasticity (Abraham and Bear 1996). Specifically, these changes would seem to limit the ability to induce further LTS. If true, this would represent a homeostatic change in the ease of inducing facilitatory plasticity, a pattern known as the BCM rule (Bienenstock et al. 1982), which serves important information-processing functions at many vertebrate synapses (for review, see Cooper and Bear 2012).

Another potential metaplasticity mechanism suggested by our data is the rapid-and-persistent up-regulation of a FMRFa transcript as well as the delayed up-regulation of a FMRFa receptor. This FMRFa transcript encodes multiple copies of the peptide neurotransmitter Phe–Met–Arg–Phe NH2 (FMRFa) (Schaefer et al. 1985). Among many functions, FMRFa has been termed a “memory suppressor” (Fioravante et al. 2006 p. 239) due to antagonistic effects on the expression of LTS memory (e.g., Abrams et al. 1984). The VC nociceptors that detect noxious shock and mediate much of the expression of LTS memory do not express FMRFa, but they do respond to FMRFa. Presumably, then, the up-regulation of FMRF transcript reflects changes in other pleural ganglia neurons, many of which are immunoreactive to FMRFa (Schaefer et al. 1985). Indeed, a number of FMRFa-expressing neurons are activated by noxious stimulation and function as interneurons in the tail and siphon-withdrawal circuits (Mackey et al. 1987; Small et al. 1992; Xu et al. 1994). Application of FMRFa to Aplysia sensory neurons antagonizes the same intracellular pathways activated by serotonin (e.g., Fioravante et al. 2006). Repeated exposure of FMRFa induces long-term depression of sensory synapses (Montarolo et al. 1988), decreases some of the same transcriptional increases produced by serotonin (Sun et al. 2001), and leads to retraction of sensory neuron processes (Schacher and Montarolo 1991). The up-regulation of FMRF signaling thus seems to represent an antagonistic or compensatory form of plasticity that helps limit the expression of LTS. It seems common for apparently unitary forms of behavioral memory to reflect a mixture of antagonistic mechanisms (e.g., Groves and Thompson 1970; Wolpaw 1997). One intriguing possibility is that the up-regulation of the FMRF system could provide an active forgetting mechanism (Berry and Davis 2014; Wixted 2004), as each further activation of the T-SWR circuit might lead to enhanced FMRFa signaling that could eventually erode the memory mechanisms induced by serotonin during training.

Materials and Methods

We report how we determined our sample size, all data exclusions (if any), all manipulations, and all measures in these studies (Simmons et al. 2012). All data for this project are posted to the Open Science Framework (https://osf.io/eaxb2); the microarray data are also posted to NCBI's Gene Expression Omnibus (Geo: GSE95596).

Animals

Animals (75–125 g) were obtained from the RSMAS National Resource for Aplysia (Miami, FL) and maintained at 16°C in one of two 90-gallon aquariums with continuously circulating artificial seawater (Instant Ocean, Aquarium Systems Inc.). Handling was as described previously (Herdegen et al. 2014b).

LTS training

A 1 d LTS training protocol (Fig. 1) was used, adapted from Wainwright et al. (2002) but with a stronger shock (90 mA versus 60 mA) and a constant-current stimulus (see Bonnick et al. 2012 for how the protocol was optimized to produce robust changes in gene expression). Training consisted of four rounds of noxious shock applied at 30-min intervals to one side of the body with a handheld electrode. Each round of shock consisted of 10 pulses (60 Hz biphasic) of 500-msec duration at a rate of 1 Hz and an amplitude of 90 mA. During the course of each shock, the stimulating electrode was slowly moved from anterior (just behind neck) to posterior (just in front of tail) and back to cover nearly the entire surface of that side of the body. Side of training was counterbalanced.

Behavioral measurement

As a behavioral outcome, we measured the duration of the tail-elicited siphon-withdrawal reflex (T-SWR) (see Walters and Erickson 1986). The reflex was evoked by applying a weak shock to one side of the tail using a handheld stimulator (60 Hz biphasic DC pulse for 500 msec at 2 mA of constant current). T-SWR behavior was measured as the duration of withdrawal from the moment of stimulation to the first sign of siphon relaxation. To characterize changes in T-SWR duration, pretest and post-test responsiveness was characterized by a series of eight responses evoked on alternating sides of the body at a 10-min ISI. Scores were split by side of stimulation (trained versus untrained) and averaged (four responses/side for each time point characterized).

The mechanisms of sensitization memory depend, in part, on the relationship between training and testing site. Site-specific sensitization occurs when the training and testing site are the same; generalized sensitization when the training and testing sites are different (Walters 1987b). These two forms of sensitization differ in the magnitude of physiological changes induced and may also be somewhat distinct in terms of their molecular mechanisms (Lewin and Walters 1999). The behavioral protocol used here has different test and training sites (tail versus side of the body) and thus measures generalized sensitization. However, we analyzed transcriptional changes in the whole pleural ganglia, which include sensory neurons mediating nociception at both the test and training sites. Thus, our transcriptional analyses represent changes evoked by both site-specific and generalized sensitization. See Herdegen et al. (2014a) for a discussion of these issues.

Isolation and processing of pleural ganglia RNA

We compared gene expression from pleural ganglia on the trained versus untrained side of the animal. Samples from two animals trained on opposite sides were pooled. This was done to control for lateralized gene expression.

To analyze transcription, pleural ganglia RNA was isolated immediately after the 24 h post-tests, a time-point representing maintenance of LTS memory. Isolation and homogenization was exactly as described in Herdegen et al. (2014b).

Reverse-transcription quantitative PCR (qPCR)

Reverse transcription was performed using Maxima cDNA kit with dsDNase (Thermo Scientific). Quantitative PCR was conducted using Maxima SYBR Green/Fluorescein qPCR Master Mix (Thermo Scientific) and the MyIQ real time PCR system (Bio-Rad). Primers were validated for correct PCR efficiency and are listed in Supplemental Table S3. qPCR samples were analyzed in duplicate and the relative amounts of each transcript were determined using the ddCT method and the Bio-Rad IQ5 gene expression analysis. All qPCR expression levels were normalized to levels of histone H4, a transcript that is stable during LTS training (Bonnick et al. 2012). For sequence confirmation, PCR products were purified and sent for sequencing to ACGT Inc.

Sample size determination

For microarray, we set a target of eight pairs of animals. This sample size was selected based on prior work, to fit the whole experiment on one slide, and to exceed the consensus recommendation of at least five biological replicates per group (Pavlidis et al. 2003; Tsai et al. 2003; Allison et al. 2006).

For qPCR validation, we set a goal of obtaining at least 10 pairs of animals. This sample size is sufficiently powered (>0.80) for large effects (d > 1) similar to what we have previously observed for regulated transcripts in this paradigm (Bonnick et al. 2012; Cyriac et al. 2013).

Quality controls

To ensure suitable samples, two quality controls were preregistered (adopted prior to data collection). First, animals had to exhibit strong learning, defined as a 1.3× change in T-SWR duration from baseline to post-test. Second, animals had to exhibit normal up-regulation of ApBiP (Kuhl et al. 1992) and ApEgr (Cyriac et al. 2013) on the trained side. We have previously used ApBiP as a quality control for 24 h analysis because it is so strongly regulated at this time point (Herdegen et al. 2014a). ApEgr was adopted as an additional control because the degree of regulation 24 h after training correlates with the degree of memory expression (Cyriac et al. 2013).

Microarray processing

We used the Aplysia Tellabs Array (ATA: GEO: GPL18666) to characterize changes in gene expression due to LTS training. This array includes 26,149 distinct probes representing all known sources of A. californica ESTs and mRNAs at the time of design (January 2012). Based on estimates from previous microarray designs (Moroz et al. 2006), the ATA should cover >50%–60% of all neurally expressed transcripts. Full details on the array design are reported in Herdegen et al. (2014b).

Microarray processing was completed by Mogene Inc. A two-color approach was used with each array hybridized to a sample from a trained or untrained animal. In half of cases, trained samples were hybridized with Cy3 and controls to Cy5; the other half we dye-swapped. Processing was exactly as described in Herdegen et al. (2014b).

Transcript annotations and de-duplication

The ATA array we used is based on ESTs. To help provide full-length sequences for the mRNAs underlying these ESTs we drew upon the transcriptome databases available at AplysiaGeneTools.org. Specifically, we used BLASTN to search for each full-length EST represented on the array in the most recent CNS transcriptome (available for download at: http://aplysiagenetools.org/wwwblast2/db/A1_CNS_merged.norm30.tgicl.fasta). Where a match could be definitively made (e < 10−10 with >20% coverage of the query) we then retrieved the full-length sequence available in the transcriptome database and used that full-length sequence for further bioinformatics analyses. Overall, 89% of the ESTs represented on the array could be definitively matched to putative full-length mRNAs. In some cases, different ESTs ended up matching the same putative full-length mRNA. We took this as evidence for duplication on the array; duplicates were removed from all downstream analyses and duplicates are not included in any of the counts reported in this manuscript.

To provide names for ESTs we matched them to the current draft of the Aplysia genome using the same approach described above (draft 3.0, annotation release 101: https://www.ncbi.nlm.nih.gov/genome/annotation_euk/Aplysia_californica/101/). Due to the nascent state of the Aplysia genome, however, only 56% of ESTs represented on the array had a definitive match with a corresponding genomic mRNA.

Throughout this manuscript, we refer to transcripts by their EST accession number and, if available, matching mRNA name. We have posted to our Open Science Framework page a full concordance between the microarray ESTs, genomic mRNAs, and AplysiaGeneTools transcripts.

Gene ontology analysis

Gene ontology analysis was conducted using standard settings with Blast2Go (Conesa et al. 2005). As the current draft of the Aplysia genome has very few functional annotations, we used tblastn to search for array probes among annotated genes from other organisms. Queries were the putative full-length mRNAs from AplysiaGeneTools or, where not available, the original EST. We used this approach to annotate all significantly up- and down-regulated transcripts. We then focused on biological process terms, identifying the most commonly occurring terms in both sets and comparing these sets for specific enrichment of terms (with an overall FDR rate of 0.05).

Statistical analyses

Behavioral responses were averaged by time point. Paired comparisons were made from baseline to post-test for each side. Standardized effect size estimates (Cohen's d) are corrected for bias (Hedges 1981) and calculated so that positive values represent an increase in response (sensitization).

Microarray data was analyzed using limma (Smyth 2005; Ritchie et al. 2015a) from the Bioconductor suite of tools (Gentleman et al. 2004) for R (Ihaka and Gentleman 1996). Our processing script for identifying differentially regulated transcripts was preregistered and is posted on the Open Science Framework. Median expression values were analyzed (Zahurak et al. 2007). These were corrected for background using the normexp+offset algorithm recommended for Agilent microarrays by Ritchie et al. (2007). An offset of 30 was selected based on inspection of MA Plots (Supplemental Fig. 2). Expression was then normalized using the loess function (Smyth and Speed 2003). Where multiple probes were used to measure the same EST or mRNA, these were averaged. Finally, trained and control expression were compared using an empirical Bayes-moderated t-test (Smyth 2004). Statistical significance was calculated using Benjamini–Hochberg correction for multiple comparisons to maintain a 5% overall false-discovery rate (Benjamini and Hochberg 1995). We used the treat function from limma (McCarthy and Smyth 2009) to conduct a stringent test for significant regulation. Specifically, rather than use a null hypothesis of no regulation, we tested for regulation statistically distinguishable from at least a 10% change in expression in either direction. We have previously found that using this type of high-stringency criterion yields very strong predictive validity in independent qPCR (Herdegen et al. 2014b; Holmes et al. 2014). Expanding beyond our preregistered analysis script we also explored the completeness of the gene list using the propTrueNull function (Ritchie et al. 2015b) and the convex decreasing densities approach developed by Langaas et al. (2005). We also explored the degree of relationship between the regulation observed 24 h after LTS training with our previous screen of regulation observed 1 h after LTS training (Herdegen et al. 2014b). We examined both the correlation between raw log-fold-change scores and the correlation once these scores are adjusted for potential measurement error using the genuine association of gene-expression profiles function (genas) in limma (Ritchie et al. 2015b)

We followed our preregistered analysis plan for quantitative PCR. A fold-change score was calculated as the ratio of expression on the trained side relative to the untrained side. For analyses, fold-change scores were log transformed (base 2). This ensures equal weight to both up- and down-regulated measures and maintains consistency with microarray analysis. Changes from control were tested using a one-sample t-test against an expected value of 0 for the null hypothesis (0 represents no change for log fold-change scores).

For ease of interpretation, fold-change scores are plotted in raw format on a log scale and are reported in text in raw format as mean fold-change (MFC) with 95% confidence intervals in brackets (Cumming and Calin-Jageman 2017). Cohen's d is reported with correction as an estimate of effect size, calculated so that positive values indicate increased expression on the trained side.