extfieldsofstudy: Biology Mathematics

s2fieldsofstudy:
Biology
Computer Science

text:
Adaptive BLASTing through the Sequence Dataspace: Theories on Protein Sequence Embedding

We theorize that phylogenetic profiles provide a quantitative method that can relate the structural and functional properties of proteins, as well as their evolutionary relationships. A key feature of phylogenetic profiles is the interoperable data format (e.g. alignment information, physiochemical information, genomic information, etc). Indeed, we have previously demonstrated Position Specific Scoring Matrices (PSSMs) are an informative M-dimension which can be scored from quantitative measure of embedded or unmodified sequence alignments. Moreover, the information obtained from these alignments is informative, even in the twilight zone of sequence similarity (<25% identity)(1-5). Although powerful, our previous embedding strategy suffered from contaminating alignments(embedded AND unmodified) and computational expense. Herein, we describe the logic and algorithmic process for a heuristic embedding strategy (Adaptive GDDA-BLAST, Ada-BLAST). Ada-BLAST on average up to ~19-fold faster and has similar sensitivity to our previous method. Further, we provide data demonstrating the benefits of embedded alignment measurements for isolating secondary structural elements and the classifying transmembrane-domain structure/function. We theorize that sequence-embedding is one of multiple ways that low-identity alignments can be measured and incorporated into high-performance PSSM-based phylogenetic profiles.

Introduction
One of the major challenges that biologists face is the ability to identify relationships between highly divergent protein sequences. Although many methods (e.g., (2;6;7)) have addressed the problem, it still remains as a challenge, as conventional sequence alignment methods often fail to obtain statistically robust measurements when sequence identity dips into the "twilight zone" (~ less than 25% identity). In general, when pairwise sequence alignments between protein sequences fall below 25% identity, statistical measurements do not provide support robust enough to identify clear phylogenetic relationships, structural features, or protein function despite intensive research in this area (2;8-11). GDDA-BLAST (Gestalt Domain Detection Algorithm -Basic Local Alignment Search Tool), originally introduced in (12), was designed to address the challenges associated with lowidentity alignments/divergence. We determined that this alignment information, when incorporated into phylogenetic profiles, is informative to our laboratory experiments at multiple scales (e.g. whole protein, single protein domain, and single amino acid) (3-5;13-16). We have used these analyses: (i) to reconstruct evolutionary histories, (ii) to identify functions in domains of unknown function, (iii) to classify structural homologues of high sequence divergence, and (iv) to inform our biochemical experimentation by isolating key amino acids important to protein function.
A phylogenetic profile of a protein is a vector, where each entry quantifies the existence of the protein in a different genome. This approach has been shown to be applicable to whole molecule (Single Profile Method), to an isolated domain (Multiple Profile Method), and to individual amino acids (17)(18)(19). GDDA-BLAST matrices are a variation of phylogenetic profiles, except in our case, a protein is a vector where each entry quantifies the existence of alignments with a PSSM(1;2). The basic idea underlying our method begins by compiling a set of PSSMs that the query sequence is compared to. These profiles can be obtained from any protein-sequence knowledge-base source (e.g. Protein Data Bank, Pfam, SMART, NCBI Conserved Domain Database (CDD) (20)(21)(22)(23)), or generated locally using PSI-BLAST (24). We employ reverse specific position BLAST (rps-BLAST (23)) to compare query and PSSMs, and have introduced multiple innovations in GDDA-BLAST. We utilize a single domain PSSM database for pairwise comparisons. Then we record and quantify all alignments between an unmodified (control), and modified query sequence. The latter is composed of two types of alignments: "seeded" and non-seeded alignments. We modify the query sequence with a "seed" from the PSSM, creating a consistent initiation site (Figure 1a-b). "Seeds" are generated from profiles by taking a portion (e.g. 10% in this study) of the PSSM sequence (e.g. from the N-terminus or C-terminus). These thresholds were used based on results from our previous studies. These "seeds" are embedded at each position of the query sequence. For example, a query sequence of 100 amino acids yields 100 distinct test sequences for each "seed". This strategy was designed to amplify and encode the alignments possible for any given query sequence. Instead of a sliding window, we utilized a sliding "seed", similar yet inverse to the embedding strategies employed by Henikoff and Henikoff(25).
Each of these modified query sequences is then aligned by rps-BLAST against the parent profile. Since BLAST algorithms are based on a "hit and extension of the hit" approach, the embedded "seed" creates consistent initiation site which allows rps-BLAST to extend an alignment even between highly divergent sequences ( Figure 1c). Next, we filter our results from rps-BLAST using thresholds of # of hits, percentage identity and percentage coverage (i.e. alignment length as a function of profile length) (Figure 1d). The output of these comparisons is a composite [product] score of either zero [when there is no significant match] or a positive value [which measures the degree of successful matching of the protein sequence to each of the profiles]. This procedure can be readily adapted to make an unbiased comparison between a series of query sequences by subjecting them to the same screening analysis with the same set of PSSM sequences as "seeds". From these results, each profile alignment above threshold is defined within the query sequence to create boundaries for our subsequent pairwise alignments. Residue scores for each amino acid position are calculated, in the absence of the embedded sequence, by scoring a value=2 for identities and value=1 for positive substitutions for all alignments above threshold. These scores are then summed from each alignment, providing a total residue score. Once scored, each query sequence (N) is represented in a vector of non-negative numbers in M dimensions (M= # of "seeds" tested) ( Figure  1e). This NxM data matrix can then be used to create a tree of relationships based on standard statistical techniques such as hierarchical clustering based on Euclidean distance or Pearson's correlation between each query sequence (Figure 1f-i).
Despite the great potential of embedded alignment strategies for answering a diverse set of biological questions, its computational cost is prohibitively expensive. This is due to the current method for generating and analyzing embedded sequences. As proteins range in length from tens of amino acids to < 8000, proteomic scale studies using GDDA-BLAST are prohibitive. To address this challenge we present here a novel sequence alignment algorithm that is as sensitive as GDDA-BLAST but orders of magnitude faster. Ada-BLAST (Adaptive GDDA-BLAST) is termed for its adaptive nature, and exploits the similarity among embedded sequences to adaptively avoid expensive computations. The algorithm begins with a modification of the query amino acid sequence at each amino acid position via the insertion of a "seed" sequence from the profile of interest. These seeds are obtained from the profile consensus sequences from NCBI's Conserved Domain Database (CDD). (C-E) Signals are collected from optimal alignments between the "seeded" sequences and profiles using rps-BLAST and are incorporated as a composite score into an N by M data matrix. (F-I) This data space can be analyzed to generate phylograms and dendrograms based on Euclidean distance measures and Pearson correlation measures of GDDA-BLAST signals, respectively.

Definitions
Let the target sequence be X and the query sequence be Y . Note that the sequence profile, to measure any query protein sequence, is called the target sequence. The length of a sequence X is denoted as | | X . Assume that | | X and | | Y are n and m, respectively. A subsequence of X from the i -th residue to the j -th residue is denoted by j i x , such that  An embedded (chimera) sequence is generated by embedding either the N-or C-terminal portion of X as a "seed", denoted by S , to every position of Y ( Figure 2). Usually % p of X (i.e., k ) is used as a seed. Thus, N-and C-terminal seeds are  C-terminal seeds from X at each position of Y ( Figure 3). Each chimera is then aligned to X using rps-BLAST. For each alignment, rps-BLAST is run independently, yielding a total of × × 2 (m # of target sequences ) BLAST executions. Note that chimera sequences are different only by the position of the seed (Figure 4). This implies that for two subsequent chimeras, much of the computation can be reused. In the Ada-BLAST approach, we re-use outputs from each step of rps-BLAST for efficient computation.
For clarity, we define the outputs of each step as follows. In the first step of rps-BLAST, we find hits between X and Y . A hit of ) , ( , where w is a word size, is  Corresponding hits between a query and a chimera sequence. This example illustrates the hits between the target sequence (X) and the query sequence (Y) can be reused for aligning a chimera sequence (C) against the target sequence (X).  (Figure 5b). If an HSP has a score high enough to trigger gapped extension, in the third step, an alignment is generated extending the HSP with gaps to both directions from a residue pair in the highest scored region of the HSP (Figure 5c-d). Note that the pair from which gapped extension is started is also referred to as a seed (24). In order to avoid confusion, we denote this as a GE starting pair to distinguish if from the embedded seed of GDDA-BLAST.

Observations of GDDA-BLAST Observation 1.
Seeding limits the search space. Since a seed provides an exact match, it is very likely that the GE starting pair is in an HSP including the seed. Also, we are only interested in the alignments including the seed because the other alignments can be found through conventional methods using the original query sequence. This limits the search space of rps-BLAST. For example, when a seed is inserted at the position 0 of a query sequence, our search space will be the region in gray (Figure 2a), and every time the seed insertion position is moved to the right, our search space is reduced as shown in Figure 2b. Note that, in the case of the chimeras with C-terminal seed, the search space is limited to the upper-left corner from the start position of the seed.

Observation 2.
Chimeras share hits. Because chimeras are the same sequence except for the position of a seed, most of their hits are conserved ( Figure 6). Therefore, we can reuse hits between X and Y to compute the alignments for any chimera sequences. Consider a chimera, C(q). Let ) (q C h be a hit obtained after the first rps-BLAST step between X and C(q). The relation between X-Y hits (i.e., h ) and X-C(q) hits (i.e., ) (q C h ) can be defined as follows: where k is the seed length, q is the seed insertion position, and w is the size of a hit. Proof is omitted because it is obvious from Figure 4. Note that we are not interested in the hits in case Figure 2a because they are out of the region of our interest by observation 1. Therefore, we only use hits in case Figure 2b for alignment.

Observation 3.
Chimeras share HSPs after ungapped extension. rps-BLAST performs ungapped extension on neighboring hits resulting in HSPs. Similar to observation 2, we can define the relationship between an HSP of X-Y (i.e., hsp ) and that of X-C(q) (i.e., ) (q C hsp ) as follows: where r is the length of the HSP. The proof is straightforward from lemma 1.

Observation 4.
Chimeras share alignment paths in gapped extension. Gapped extension in rps-BLAST starts at a GE starting pair that is a central residue pair in the highest scoring segment of any HSP whose score is sufficiently high. Different alignments can be generated if gapped extension is performed from different GE starting pairs and there is no guarantee that the same GE starting pair is selected for different chimera sequences. However, as shown in Figure 4 (boxed residues represent the seed), if a portion of a target sequence is conserved in a query sequence, it is very likely that the conserved region is aligned for multiple neighboring chimera sequences. We exploit this property to speed up the alignment process.

Observation 5.
Not every chimera produces a useful alignment. Even though a seed provides artificial matches, it cannot be extended if there are insufficient HSPs nearby to connect to. Therefore, we can significantly reduce the computational complexity of the alignment process by inserting seeds only to a limited number of query positions that are likely to be extended. Hence, in Ada-BLAST, we align the query and the target sequence first and then compute the seed insertion positions from the alignment result before aligning the chimera sequences.

Ada-BLAST Algorithm
Ada-BLAST works in four basic steps as shown in Figure 5. First, we find conserved regions by generating non-overlapping local alignments between the query and the target sequence (26). We call these partial alignments. Second, for each partial alignment from step 1, seed insertion positions are determined. Third, we produce final alignments including the seeds. Finally, we filter out nonsignificant alignments using quality measures such as the % coverage and % identity of the alignment to the corresponding PSSM.
Step 1. Find multiple non-overlapping local alignments. Haung et al (26) proposed an algorithm to find multiple non-overlapping local alignments between two sequences. We adopted this algorithm to generate partial alignments between the query and the target sequence. For any target sequence, we can use any given substitution matrix (e.g. BLOSUM62, BLOSUM45, PAM30, etc) in our algorithm.
The local alignments are found as follows. First, hits between the query and the target sequence are found. A hit is three consecutive residues with score larger than a threshold (i.e. minimum word score). For each hit, to generate an HSP, ungapped extension is performed until the score drops below a threshold (HSP drop-off score). When a new HSP is constructed, hits involved in the HSP are removed to prevent subsequent HSPs extending over them. This ensures that all local alignments produced later are non-overlapping. We keep only the HSPs with score larger than a threshold (minimum HSP score). Gapped extension is then performed for each HSP to generate partial alignments.
We keep only the partial alignments whose lengths are greater than threshold (i.e. minimum partial alignment length). Minimum partial alignment length is determined proportionally to the length of the target sequence. Since a partial alignment is a locally best alignment, it is not likely that a seed will be extended further than the end position of the partial alignment. In the final step, a final alignment will be filtered by the coverage of the alignment over a target sequence. This pre-filtering on partial alignments can remove seed insertion positions where a seed cannot be extended to final alignments with sufficient coverage.
Step 2. Select seed insertion positions. In this step, seed insertion positions in the query sequence are selected given the partial alignments obtained from step 1. As discussed in observation 1, a final alignment is generated by extending a seed from its end positions. Since a partial alignment is a locally optimal alignment, the extension of a seed can be converged to the partial alignment if a seed is inserted nearby and the score of the path is high enough. Because of the relatively high penalty of gaps used in sequence alignment methods, an alignment is usually generated with HSPs connected with small numbers of gaps in between. In addition, the score of partial alignment on either side of the gaps must be higher than the gap penalty [17]. The gapped extension usually starts from the seed to the partial alignment because the score of the alignment with the seed is typically much higher than that of the partial alignment. For this reason, we can compute the seed insertion positions simply with the score of a seed and the distance from the seed to the partial alignment.
Given a seed S of score Score(S), the maximum gap G(S) (i.e., the distance from the seed to the partial alignment) is computed as follows: , where GOP and GEP are gap opening penalty and gap extension penalty, respectively. Given a query sequence Y and a partial alignment ( b a x , , d c y , ), the query position q is subject to embedding a seed of length k as follows: (1) for N-terminal seed: , and (2) for C- , no C-terminal seed is inserted. The idea of maximum gap has been described previously for connecting HSPs with gaps(9).
Step 3. Generate final alignments with a seed. For each query position q identified in step 2, we perform an alignment with the seed S inserted in the respective position to generate the final alignments. The final alignments are generated by running dynamic programming starting at the end position of the seed, ( | S | -1, q), and proceeding to ( | X | -1, | Y | -1). Since we are working with highly divergent sequences producing low-identity alignments, it is reasonable to consider the scenario that a longer alignment with lower score can be biologically more meaningful than a shorter alignment with higher score (11). Motivated by this observation, during the alignment, we adjust an alignment score with respect to the length of the alignment as follows: is a cell in the dynamic programming matrix, and Step 4. Filter out non-significant alignments. Not all alignments produced from the previous step are informative. In this step, we prune out insignificant alignments using the metrics, % coverage and pairwise identity. Given an alignment ( b a x , , d c y , ), the % coverage of the alignment to a target sequence is calculated as follows: where X is a target sequence. The pairwise identity considered here is the identity of the alignment excluding the matches in a seed: i.e. pairwise identity= 100 # × α matches of , where α = # of matches + # of mismatches + # of gaps in an alignment excluding the seed. If the pairwise identity and % coverage of a final alignment are greater than the thresholds, minimum identity and minimum coverage, the alignment is returned to the user.

Preparation of function or structure-specific PSSM set
To generate the PSSM set for a specific protein function or structure fold, we first collected protein sequences which are known to be related to the function or structure of our interest. For the PSSM set, we generated PSSMs with the collected sequences or the sequences expanded from the collected sequences using PSI-BLAST (24). For expansion, each collected sequence is searched against NCBI NR database by PSI-BLAST (with the option of -e 1X10 -3 -h 1X10 -6 ). Among the returned sequences, we filtered out the sequences whose pairwise identity to the query sequence is more than 90% and redundant sequences in the set. And, for PSSM generation for those expanded sequences, PSI-BLAST (with the option of -h 1X10 -6 ) was used again. All PSSM sets which have been used for our test will be provided upon request.

Experimental setup and data sets
Both GDDA-BLAST and Ada-BLAST were implemented in C, and compiled for both UNIX and Windows environments. GDDA-BLAST utilizes rps-BLAST in NCBI BLAST 2.2.15 package to compute the alignments. In order to validate the approach, we tested both for the execution time and the accuracy. The execution time experiment was conducted in a dedicated machine with 1.8GHz Intel Core TM 2 duo processor and 2GB memory running Windows Vista. The experiment for accuracy was performed on a server with eight Dual-core 2.4 GHz AMD Opteron processors and total of 32G memory running Linux. Note that for the execution time experiment we used a less-equipped dedicated machine instead of the server shared by others in order to measure the execution time more accurately.

Execution time and accuracy
To compare the execution time of GDDA-BLAST and Ada-BLAST , we ran both methods with 602 query sequences randomly chosen from the SABmark twilight zone set (27) and 51 target sequences randomly selected from the CDD database. Figure 7 shows the per-query execution time when a given query is run against the 51 PSSMs in the library. The lengths of the 602 query sequences range from 34 to 759 amino acids. Note that the running time of GDDA-BLAST increases linearly as a function of query sequence length. Conversely, Ada-BLAST shows much better scalability with respect to the size of query because it inserts a seed only at the positions where the seed is likely to be extended. Moreover, the performance gain is maximized when the two sequences compared are of low identity because the number of likely seed-insertion positions is limited. This makes Ada-BLAST an attractive alternative for the alignment of highly divergent sequences. Overall, Ada-BLAST is 19.3 times ( ± 15.29 S.D.) faster than GDDA-BLAST on average while it achieves more than 100 times speed-up in many occasions.
We performed Receiver Operating Characteristic (ROC) curve analysis [4]  In this evaluation, we used 534 sequences out of 61 sequence groups randomly selected from the SABmark Twilight zone set, as described previously(1). Pearson's correlation is used for GDDA-BLAST and Ada-BLAST and sequences with K highest Pearson's correlations to the query are returned. The difference in accuracy of GDDA-BLAST when compared with Ada-BLAST is negligible.

Characterization of Secondary Structural Elements
In order to determine the information content contained in a pure population of embedded alignments we analyzed the structurally resolved (X-ray Crystallography) transmembrane protein, Bovine Rhodopsin (PDB: 1F88) (28). Figure 8A depicts the output of rps-BLAST (e-value threshold 0.01) for the domain architecture of 1F88. Notably, rps-BLAST returns overlapping alignments for 5 different PSSMs defined as Serpentine type 7TM domains. Our theories on structurally/functionally related PSSM libraries predict that additional information below accepted statistical thresholds can be utilized to define, with higher resolution, domain boundaries and secondary structural elements. This hypothesis was tested and the performance was evaluated against rps-BLAST and Hidden Markov Models (TMHMM) (29). For this experiment we considered multiple running conditions for rps-BLAST: (i) run rps-BLAST against the CDD database, and among the alignments reported, consider only those PSSMs which we have previously annotated as transmembrane by keyword (CDD 539), (ii) run rps-BLAST against a database of PSSMs (n= 24,378) derived from expanding all of the sequences contained in the original 539 PSSMs using PSI-BLAST (integral lipidbinding database, ILB DB) (see Methods), and (iii) run rps-BLAST against these databases and slide the e-value threshold to less statistically significant levels (0.01, 100, or 1xe 10 ).
The graphs in Figure  8b show the alignment results which were smoothed by adjacent averaging with a 7 amino acid window. For each case the transmembrane probability determined by TMHMM is shown on the left axis (magenta). The right axis represents a positional score for Ada-BLAST and rps-BLAST conditions. The positional score was quantified in the following manner. For each positive PSSM, the alignment boundaries as determined by the overlapping alignments obtained from either Ada-BLAST or rps-BLAST. These regions were extracted and realigned by the Smith-Waterman algorithm with a BLOSUM62 and BLOSUM45 substitution matrix. Identical residues were scored value=2 and those which were positive (non-identical but conserved) were scored value=1. This process was repeated for all positive PSSMs and the results summed for each amino acid in the protein. Currently, we do not know what substitution matrix is the most optimal so we averaged the positional scores derived from BLOSUM62 and BLOSUM45. The positional results in Figure 8B (i-v) were normalized to zero by subtracting the average positional score across the protein length from each point and then each amino acid position was subjected to adjacent averaging. These results demonstrate that: (i) using the expanded ILB DB increases the signal-to-noise ratio using either rps-BLAST or Ada-BLAST, (ii) Ada-BLAST has a ~10-fold increase in signal compared to the largest results obtained from rps-BLAST (1xe 10 ), and (iii) the positional data, even at the highest statistical limit of rps-BLAST tested with the expanded ILB DB, accord with results obtained from TMHMM. . Compared with rps-BLAST, Ada-BLAST is more refined with respect to the annotation of alpha-helices. Moreover, this data demonstrates that less statistically valid alignments (e.g. e-value 0.01 vs 100) are still informative for detecting domain boundaries and outperform lower thresholds. The full-length structure of Rhodopsin is shown (dimer) as well as an inset of the C-terminus which is composed of three small helices with the last one folding parallel with the membrane (it itself is not believed to be transmembrane).
It is reasonable to consider that amino-acids within transmembrane spanning helices will be more conserved than the intervening loop residues. Support for this hypothesis is presented in Figure 9, wherein we report our results when compared to the known structural elements of 1F88 as obtained from X-ray crystallography (28). The full-length structure of Rhodopsin is shown in the bottom right (dimer) as well as an inset of the C-terminus. For this analysis, we used Origin Lab 7.5© to perform smoothing (Fourier-transform point=8) and discontinuous baselining of the positional data. This correction was performed by baselining to every local minimum across the entire curve. The structural features are annotated with droplines (Cyan= Beta pleated sheets, Green=helices, loops not shown). In all cases, the curves obtained from the positional data can be correlated to the known structural elements.
While none of the prediction methods accurately model all of the crystal structure, we observe several interesting features. For example, several of the membrane-spanning helices are interrupted by loop regions that are not identified by TMHMM. Indeed, the C-terminus of 1F88 contains 3 small helices, the last of which is a bent-helix believed to be parallel to the membrane (aa 288-348, Figure  9 inset). Both rps-BLAST and Ada-BLAST detect these smaller helices with Ada-BLAST having the highest signal. Another region of interest is contained between aa 91-111, which is a loop in the crystal structure, but is predicted to be a short helix by rps-BLAST and Ada-BLAST. We theorize that this loop may be, under native conditions, a bent-helix similar to other regions in the protein.

Classification based on Hierarchical Clustering
We hypothesize that the alignment information described above can be used to classify protein structure/function when encoded into a phylogenetic profile. In general, similarity measurements between two protein sequences are done by directly aligning the two sequences one against the other. However, GDDA-BLAST and Ada-BLAST do not compare the two query sequences directly. Instead, we compute phylogenetic profiles for each query sequence which provides a quantitative platform for sequence comparison(1;2;4). In addition, we have made several modifications to rps-BLAST, such that its data can be encoded into phylogenetic profiles.
A phylogenetic profile of a query is a vector of length equal to the number of target sequences such as ) , , , where m is the number of target sequences. Each dimension of the vector, j b , represents the alignment score for the query sequence against the th j target sequence and computed as where j h is the number of alignments between the query sequence and the target sequence after filtering (in step 4), j s is the average %identity of the alignments, and j v is the maximum coverage among the alignments. Note that j h can be up to the length of the query sequence (i.e., #of embedded seeds) because we align the embedded sequence of the query against the target. In the case of rps-BLAST, the calculations are the same except for the presence of embedded sequence.
We use Pearson's correlation coefficient to measure the similarity between two query profiles. Pearson's correlation coefficient between two profiles A and B is computed as follows (30):  Figure 10 Classification based on Heirarchical Clustering. 74 sequences representing multiple classes of transmembrane containing proteins were hierarchically clustered by the microarray analysis algorithm Cluster (settings: Complete Linkage and Correlation (centered) similarity metric) and visualized utilizing the Treeview algorithm. Red lines represent the correlation scores derived from the analysis. Alignments for the 24,378 ILB DB PSSMs were derived by either (a) Ada-BLAST, or (b) rps-BLAST.

Fig 10a
To test the efficacy of phylogenetic profiles built using either a pure population of embedded alignments or alignments over a range of rps-BLAST e-value thresholds, we curated a set of transmembrane containing proteins from a range of different protein families (voltage-gated Ca 2+ , K + , and Clchannels, calcium-activated K + channels, cyclic-nucleotide gated channels, transient receptor potential channels (TRPs), receptor tyrosine kinases, G-protein coupled receptors (GPCRs),  TRPC   TRPM   TRPV2   TRPV   TRPA1  TRPM8   TRPP   TRPM5  TRPM2   TRPML   TRPM4   Fig 10b transporters, and exchangers). Each of these 74 sequences was analyzed using the expanded ILB DB with Ada-BLAST and rps-BLAST. The results from this analysis were encoded into a NXM matrix using the scoring scheme ([# of hits] x [mean % identity] x [max % coverage]) and then subjected to hierarchical clustering. In Figure 10 a-b, we report the conditions that achieve the highest degree of classification (Ada-BLAST and rps-BLAST e= 0.01 respectively).
While neither of these classifications are perfect, both are robust for pairing related sequences; rps-BLAST has the highest Pearson's correlation, while Ada-BLAST has a more robust topology. Specifically, Ada-BLAST classifies the TRP channels, GPCRs, and transporters/exchangers into their correct group, while rps-BLAST does not. In data not shown, we observe that both the Pearson's correlation values and the overall topology is compromised when measured at e= 100 and e= 1xe 10 in rps-BLAST. Taken together, these data demonstrate that a pure population of embedded alignments can differentiate between proteins having similar structure and diverse function given the appropriate PSSM database.

Discussion
In this study, we provide evidence for our theories that: (i) sequence embedding amplifies lowidentity alignments and that these alignments are distinct from those derived by simply scaling evalue thresholds, (ii) low-identity alignments contain information which can inform protein structure/function annotation before and after encoding into a phylogenetic profile, and (iii) PSSM libraries which are constructed from key-word searches of the CDD database, expanded using PSI-BLAST and then made into a specific database can be used to classify proteins based on their structure/function. There are several implications which can be drawn from these findings.
The biochemical characterization of membrane spanning proteins is challenging. Likewise, structural studies of this protein class can be fraught with artifacts introduced by crystallization and/or lack of appropriate co-factors (e.g. lipids). Our results demonstrate that signals derived by Ada-BLAST can provide structural information that is distinct from Markov Models of transmembrane regions. It is tantalizing to consider that these measurements may isolate potential discrepancies in protein structure due to non-native conditions. Importantly, when encoded into a phylogenetic profile, these same data can be utilized to functionally classify transmembrane containing proteins.
The power underlying protein phylogenetic profiles is centered on the construction of PSSMs representing specific folds and/or activities. Indeed, results from this manuscript and also our companion study (Ko et al, Physics Archives November 2009) demonstrate that PSSMs libraries generated using specific folds can accurately identify homologous folds. Moreover, PSSM libraries generated for a specific activity can accurately identify homologous functions in proteins of diverse structure, as well as differentiating activity within a specific fold. Within the present study we generated transmembrane PSSMs curated by key-word searches; however, it is likely that additional refinement of this set (e.g. to remove non transmembrane regions from PSSM) will lead to an increase in the signal/noise ratio, better annotation of transmembrane spanning domains, secondary structure prediction, and more robust classification. We are actively pursuing our working hypothesis that this approach will work for any class of protein domains.
Our results support the idea that statistical thresholds are often too stringent in domain detection algorithms. For example, rps-BLAST does not report a channel domain alignment in human TRPV channel (gi|22547180) at statistical limits. In the present study we found that additional information contained in alignments well below accepted statistical thresholds can be utilized to inform domain boundaries and secondary structural elements. However, this was not the case when these same data were hierarchically clustered. Future analysis on sufficiently large datasets is required to identify and optimize the multiple variables which can identify highly-divergent yet informative alignments. Nevertheless, we propose that there is a wealth of information below statistical values which can aid researchers in annotating protein structure/function.
In conclusion, we propose that future work aimed at (i) creating comprehensive and refined PSSM libraries and (ii) exploration of sequence embedding at the level of the PSSM (COBBLER, (25)) and within the query (Ada-BLAST), will exponentially increase the functional annotation of all classes of proteins across taxa. Such an advance would have broad impacts on human health and disease, as well as basic science endeavors. Indeed, since Ada-BLAST performs in the "twilight zone" of sequence similarity, this approach can be harnessed to decode the most challenging protein datasets, and scaled up to screen proteomes and the vast quantities of sequences being obtained from metagenomic studies. Outside of biological questions, the theories behind these algorithms are likely to have applications in many other fields that use pattern-based prediction algorithms.


Genomic instability in quartz dust exposed rat lungs: Is inflammation responsible?

Exposure to quartz dusts has been associated with lung cancer and fibrosis. Although the responsible mechanisms are not completely understood, progressive inflammation with associated induction of persistent oxidative stress has been discussed as a key event for these diseases. Previously we have evaluated the kinetics of pulmonary inflammation in the rat model following a single intratracheal instillation of 2mg DQ12 quartz, either in its native form or upon its surface modification with polyvinylpyridine-N-oxide or aluminium lactate. This model has been applied now to evaluate the role of inflammation in the kinetics of induction of DNA damage and response at 3, 7, 28, and 90 days after treatment. Bronchoalveolar lavage (BAL) cell counts and differentials as well as BAL fluid myeloperoxidase activity were used as markers of inflammation. Whole lung homogenate was investigated to determine the induction of the oxidative and pre-mutagenic DNA lesion 8-hydroxy-2-deoxy-guanosine (8-OHdG) by HPLC/ECD, while mRNA and protein expression of oxidative stress and DNA damage response genes including hemeoxygenase-1 (HO-1) and apurinic/apyrimidinic endonuclease (APE/Ref-1) were evaluated using Western blotting and real time PCR. Isolated lung epithelial cells from the treated rats were used for DNA strand breakage analysis using the alkaline comet assay as well as for micronucleus scoring in May-Gruenwald-Giemsa stained cytospin preparations. In the rats that were treated with quartz, no increased 8-OHdG levels were observed, despite the presence of a marked and persistent inflammation. However, DNA strand breakage in the lung epithelial cells of the quartz treated rats was significantly enhanced at 3 days, but not at 28 days. Moreover, significantly enhanced micronucleus frequencies were observed for all four time points investigated. In the animals that were treated with the PVNO modified quartz, micronuclei scores did not differ from controls, while in those treated with the aluminium coated quartz intermediate effects were found. These findings were in line with the kinetics of inflammation and epithelial proliferation in the rat lungs for the different treatments. Notably, a highly significant correlation was observed between neutrophil numbers and micronucleus frequencies, indicative for a role of inflammation in eliciting genomic instability in target cells of quartz-induced carcinogenesis. Our ongoing investigations focus on the evaluation of the causality between both in relation to quartz exposure.

Introduction
Inflammation has been recognised as an important factor in cancer development and progression [1,2,3]. Chronic inflammation is also considered to be crucial in driving lung cancer formation after inhalation of quartz as well as poorly-soluble low toxicity dusts. Herein, specifically a mechanism is considered to be involved whereby inflammatory cells release excessive amounts of reactive oxygen species (ROS) which are known to be able to cause oxidative DNA damage and mutagenesis [reviewed in 4]. This mechanism can be referred to a secondary genotoxicity, in contrast to primary genotoxicity, which is considered to involve a direct genotoxic action of a particle [5]. Previous studies have demonstrated the carcinogenicity of quartz as a variably entity [6]. We and others have shown that especially the surface reactivity of the particles contributes to their biological response, since surface modifications by polyvinylpyridine-N-oxide (PVNO) or aluminium lactate (AL) are able to reduce the inflammatory and genotoxic potential of DQ12 quartz particles [7,8,9,10]. The present study was undertaken to evaluate the relation between the kinetics of quartz-induced lung inflammation and genotoxicity in vivo at 3, 7, 28 and 90 days after single intratracheal (i.t.) instillation. Inflammation was characterised by inflammatory cell counts of the bronchoalveolar lavage (BAL) and differentials. For the evaluation of genotoxicity three independent methods were used. In whole lung tissue homogenate we determined the formation of the oxidative and pre-mutagenic lesion 8-hydroxydeoxyguanosine [11]. Two further measurements of genotoxicity were determined specifically in the lung epithelial cells, isolated from the treated animals. These cells were specifically analysed because of their relevance as target cells in pulmonary carcinogenesis [8]. For the determination of DNA strand breakage the alkaline comet assay was used. This assay is a sensitive genotoxicity indicator assay that allows for quantification of DNA-damage in individual cells [12]. The isolated lung epithelial cells were also evaluated for micronuclei (MN), which are small extranuclear bodies, arising in dividing cells from chromosome fragments or a whole lagging chromosome. MN analysis in isolated lung epithelial cells from rats has been successfully applied for the evaluation of in vivo genotoxicity of hard metal particles and carbon nanotubes after i.t. instillation [13,14].

Materials and methods
2.1. Particle characteristics and coating procedure Dörentruper quartz (DQ12) was used as native quartz as well as for the preparation of surfacemodified quartz samples. The procedure for particle coating with aluminium or PVNO and quality control evaluation (i.e. determination of coating efficiency and size distribution) have been described elsewhere [8,9,10].

Animals and treatment
Female Wistar rats (Janvier, Le Genest St. Isle, France) were maintained in an accredited on-site testing facility, responding to the guidelines of the Society for Laboratory Animals Science (GV-SOLAS). At a age of eight weeks (220±1g) animals were i.t. instilled following anesthetization (Isofluran, Essex Pharma GmbH, Munich, Germany) using 400µl volume containing phosphate buffered saline (PBS) as vehicle control, 2mg DQ12 quartz or DQ12 (2mg) coated with either PVNO or with AL. At days 3, 7, 28 and 90 post-instillation, animals were sacrificed by deep anesthetization with pentobarbital (50 mg/kg body weight) followed by exsanguination via the A. abdominalis.

Broncho-alveolar lavage and analysis
Lungs of 5 animals per treatment group and time point were lavaged as described previously . The BAL fluid was spun (500g, 10min, 4°C), and cells were collected for cell counting and May-Grunwald-Giemsa (MGG)-stained cytospin preparations. Inflammatory and cytotoxic data are published elsewhere [9].
2.4. DNA isolation and analysis of 8-hydroxy-2'-deoxyguanosine by HPLC/ECD Lung tissue was removed, chopped, aliquots were snap frozen in liquid nitrogen and stored at -80°C until later measurement of 8-OHdG using high performance liquid chromatography with electrochemical detection (HPLC-ECD) [15]. Values are expressed as the total amount of 8-OHdG per 106 nucleotides.

Lung epithelial cell isolation
In order to investigate the specific target cell for particle-induced lung tumours, lung epithelial cells were isolated according to the method of Richard and colleagues [16] as modified by Knaapen and coworker [8]. A part of the obtained epithelial cell fraction was used to prepare MGG-stained cytospin preparations in order investigate the appearance of MN. Another part was used for the evaluation of DNA strand breakage by the alkaline comet assay.
The purity of the isolated cell fraction, as determined by microscopy evaluation of the MGGstained slides was found to be 63 -75 % epithelial cells. Other cell types were neutrophils (4 -18 %) and macrophages (13 -20 %).
2.6. Comet assay DNA strand breakage was investigated immediately after cell isolation using the alkaline comet assay as described previously [8] according to the guidelines published by Tice and Colleagues [17]. For each animal three independent slides were prepared, and on each single slide 50 cells were randomly evaluated. Cells were categorised according to presence or absence of a tail. The comet assay analysis was performed for 2 time points, i.e. at 3 days and 28 days.

Micronucleus evaluation in lung epithelial cells
In order to determine the proportion of the micronucleated cells, in each MGG-stained slide 1000 lung cells were counted at x1000 magnification with immersion oil using a Zeiss Axioscope 2 microscope. The identification and scoring of the MN was done according to Tolbert and colleagues [18]. Characteristics of each evaluated cell are an intact cytoplasm and an intact nucleus with a smooth and distinct nuclear perimeter. A cellular structure was then identified as a MN if its size was less than a third of the diameter of the associated nucleus, but large enough to discern the round shape and a comparable staining intensity to the nucleus.

Statistical analysis
Data are expressed as mean ± SD. Statistical analysis was performed using SPSS version 15.0 for Windows, post hoc Tukey-HSD. Correlation between MN formation and percentage of PMN was tested using Person's r-test. Differences were considered as statistically significant when p < 0.05.

Results
Investigation of whole lung tissue homogenate showed no enhanced 8-OHdG/dG ratios in the animals after treatment with native quartz at all investigated time points (Figure 1). Surprisingly, at all time points 8-OHdG levels tended to be higher in the lungs from rats that were exposed to the coated quartz preparations than from those that were treated with the original quartz sample. At the 7 day time point, 8-OHdG levels were significantly enhanced in animals that were treated with the PVNO-coated quartz sample (p<0.05, Fig. 1). The comet assay was used for the determination of DNA strand breakage in freshly isolated epithelial cells. The measurements were performed at the 3 day and the 28 day time points (Table 1). Table 1: Percent of isolated lung epithelial cells showing DNA damage as determined by the alkaline comet assay. Cells were isolated from the lungs of rats that were exposed to a single dose of 2 mg noncoated DQ12 or DQ12 coated with PVNO or AL (* p< 0.05). Three days after quartz instillation a significant higher number of cells with DNA damage were observed, when compared to the control animals. Particle surface modification with PVNO and AL inhibited the ability of quartz to induce these DNA damage. After 28 days, these differences were not present. In fact, the DQ12 treated animals tended to show the lowest DNA damage at this time point.
As a further marker for the genotoxic potential of the quartz particles, the appearance of MN were evaluated at MGG-stained cytospin preparations of lung epithelial cells. In quartz-treated animals significantly higher MN % were observed at all times points ( Figure 2). Interestingly, for animals that had received AL-coated quartz also enhanced MN frequencies were found at the 28 and 90 day time points. In contrast, animals that were treated with the PVNO-coated quartz MN % in the lung epithelial cells isolates were not increased. In order to evaluate the relation between the inflammation as induced by the different quartz preparations and genotoxicity, correlation analyses were performed between inflammatory cell counts from bronchoalveolar lavage and the MN frequencies as determined in the lung epithelial cells. For each of the four investigated time points there was a significant linear correlation between the % of neutrophils and the MN frequencies when data were considered on the single animal level (3,7,28 and 90 days, all: p<0.01). The relation between neutrophilic inflammation and epithelial lung cell genotoxicity was also observed on the group level for the entire study ( Figure 3). This effect was also clearly exhibited on the single animal level (n = 76; R 2 = 0.46; p < 0.0001).

Discussion and conclusion
The major goal of our study was to evaluate the relation between inflammation and genotoxicity, considered to be a crucial mechanism of particle-induced carcinogenesis. Therefore, we evaluated the kinetics of lung inflammation after single i.t. instillation of different preparations of quartz in rat lungs in relation to various genotoxicity readouts. The application of two different surface-modifications of DQ12 quartz (i.e. PVNO and aluminium) along with the original DQ12 sample resulted in levels of inflammation which largely differed in extent as well as over time [9]. In this study we have now examined over the time different markers of particle-induced genotoxicity. Being the most-well investigated oxidative DNA lesion, the induction of 8-OHdG was measured in whole lung tissue homogenates by HPLC/ECD [11,15]. Increased 8-OHdG immunoreactivity in lung alveolar cells has been previously shown by Seiler and colleagues [19] three weeks after a single i.t. instillation of 1.5 mg quartz in rats, while 3 days after instillation such an increase was absent. Surprisingly, no significant 8-OHdG induction was observed by the original quartz sample in the current study despite a persistent inflammation, neither with HPLC/ECD nor with supportive investigations by immunohistochemistry (data not shown). The divergence between the different in vivo quartz studies may have resulted from use of different doses, exposure times or quartz samples, i.e. as an effect of a "batch to batch variation" [6]. Also the possibility of an artificial induction of 8-OHdG during DNA extraction and processing as required for HPLC/ECD analysis, leading to high background levels would neither explain positive findings by others [e.g. 11] as well as our unexpected effects in the animals which received PVNO-coated quartz samples. In contrast to the lack of detectable increases in 8-OHdG upon DQ12 treatment, a significant increase in DNA strand breakage of lung epithelial cells from quartz-treated animals could be measured three days after quartz instillation. This difference might be explained by the specific isolation of the lung epithelial cells and the high sensitivity of the comet assay method. Importantly, DNA-damage was not observed in the animals that were treated with the surface modified quartz samples. As these data were related to the extent of neutrophilic inflammation as induced by the different sample treatments for this time point, they are in support of the causal role of inflammation in particle-induced genotoxicity [8]. However, 28 days after the instillation, no difference in DNA strand breakage could be detected in the lung epithelial cells from the different treatment groups. These findings may be explained by an induction of DNA repair, as indicated from our previous investigations were we observed enhanced expression of the oxidative DNA damage repair protein apurinic/apyrimidinic endonuclease/redox factor 1 (APE/Ref1) [10]. However, the role of epithelial cell proliferation and differentiation known to occur after quartz treatment should also be considered. In this regard, the third genotoxicity assay evaluated in our current study, i.e. MN in the lung epithelial cells, are considered to be the most relevant: While the comet assay is a measure of the bare presence of a potentially broad spectrum of DNA lesions which may -or may not -lead to a mutation, a MN represents the product of a clastogenic or aneugenic effect and cell proliferation. The relevance of MN test has been further supported by recent observations that demonstrated the predictive value of their increases in peripheral blood lymphocytes for lung cancer risk in humans [20].
In the present study, MN frequencies were found to be enhanced in lung epithelial cells up to the period of 90 days after quartz instillation. In the knowledge that the quartz treatment in our study has been shown to lead to a persistent increasing inflammation over this time period [9], it may be proposed that the observed in vivo genotoxicity results from (a) yet to be identified mechanism(s) operating in inflamed lung tissues. This is further supported by observed gradients in inflammation and MN in the lungs of the rats that received AL-and PVNO-modified DQ12 respectively, visualised by a highly significant correlation between the percentages of neutrophils and MN in the over all study. Although previous in vitro studies by us and others have demonstrated that both neutrophils and quartz particles can cause genotoxic insult to epithelial cells [reviewed in 5] our present data are in further support of the secondary mechanism of genotoxicity. Noteworthy in this regard, the MN frequency in epithelial cells of the control rats as well as those that were treated with PVNO-coated DQ12 was about 0.5 -1 % indicative of the absence of an increased genotoxicity up to an approximate 30 % of neutrophils in the lung as observed with the DQ12-PVNO group.
In conclusion, our study demonstrates that modification of the reactive surface of quartz particles can have significant impact on pulmonary inflammation and associated genotoxicity and therefore provides further clues for the observed variable cancer risk of quartz exposure [6]. Moreover, our data are the first to show a marked in vivo association between neutrophilic inflammation in the lung and the formation of MN in lung target cells for (particle-induced) lung carcinogenesis. Additional experiments are needed to establish the causality of this relation.


Carbon Nanotube Membranes for use in the Transdermal Treatment of Nicotine Addiction and Opioid Withdrawal Symptoms.

Transdermal systems are attractive methods of drug administration specifically when treating patients for drug addiction. Current systems however are deficient in therapies that allow variable flux values of drug, such as nicotine for smoking cessation or complex dosing regimens using clonidine when treating opioid withdrawal symptoms. Through the use of functionalized carbon nanotube (CNT) membranes, drug delivery to the skin can be controlled by applying a small electrical bias to create a programmable drug delivery system. Clearly, a transdermal patch system that can be tailored to an individual’s needs will increase patient compliance as well as provide much more efficient therapy. The purpose of this paper is to discuss the applicability of using carbon nanotube membranes in transdermal systems for treatment of drug abuse.

Introduction
Drug abuse, whether alcohol, tobacco, prescription or illicit, continues to burden society, economics, medical care systems, and the family dynamic. The World Health Organization (WHO) estimates that globally 185 million people are illicit drug users and an astounding 1.3 billion individuals are smokers. 1 Smoking related deaths, estimated at nearly 5 million worldwide, continue to be the leading preventable cause of death. Additionally, in the United States alone, $167 billion in annual health-related economic loss can be attributed to smoking. 2 Although there are roughly only 200,000 deaths reported to be directly associated with illicit drug use, the detriment to family and society, and the economic costs associated with crime, law enforcement, medical care, and drug treatment and rehabilitation, are immense. 3,4 While current treatments for cigarette smoking and opioid addiction do aid in successful abstinence from the individual's vice, there is great need and room for improvement. The best of smoking cessation therapies reports a mere 19% success after one year, while only 50% of individuals who enter clinics for opioid addiction remain for the duration of the treatment, of which 65% return to using within 4 weeks. 3,5 Clearly, a therapy that better controls cravings or abates withdrawal symptoms to decrease the likelihood of relapse is needed. Carbon nanotube (CNT) membranes present the opportunity to create a transdermal patch that can vary its rate of delivery throughout its application to the skin to attain therapeutic plasma levels and plasma profi les of a specifi c drug. The purpose of this report is to discuss the use of carbon nanotube membranes in a transdermal system to aid in abstinence from cigarettes and withdrawal from opioids.

Nicotine and Clonidine as Therapies
Conventional nicotine replacement therapy (NRT) has been on the market since the early 1990s. Continuous research and development of the therapy has expanded the delivery from early reservoir transdermal patches to matrix systems, gums, inhalers, lozenges, and nasal sprays. 6 As the name implies, NRT works to replace nicotine, the agent in cigarettes responsible for dependence, with a pure pharmaceutical source, thereby eliminating the toxic chemicals found in tobacco smoke from entering the body. Furthermore, NRT prevents withdrawal symptoms which induce the desire to smoke, thereby increasing a smoker's ability to put up with the environmental and social cues associated with smoking. Though current NRTs do aid in smoking cessation when compared to individuals who try to quit using non-pharmacological treatment or are given placebos, the success rate is still an unsatisfactory 17%. 3 A representative daily nicotine plasma profi le of a heavy smoker can be seen in Figure 1a. It consists of a series of peaks and troughs spurred by the smoking of each cigarette with rapid, less than fi ve minutes, attainment of the maximum nicotine level. Slight accumulation can be seen over the course of the day, but the high spike of nicotine is always obtained with the smoking of a new cigarette. 7 NRT's low success rate is primarily due to the fact that no current system can actually mimic the rapid attainment of high nicotine plasma levels associated with smoking a cigarette. 3,8 Nasal spray offers the fastest time of delivery to its maximum concentration, but does not reach the level required to suppress a nicotine craving of a smoker, and like the gum, lozenge or inhaler it requires multiple doses throughout the day which often leads to issues with patient compliance. The current transdermal patch is a once daily application, however it cannot safely deliver the craving sustaining dose, nor can the system deliver it in a rapid manner, mimicking the smoking plasma profi le. 6 Clearly, a delivery system that can rapidly and safely deliver nicotine at a level equal to that of a cigarette, followed by a return to baseline levels, and then later be reactivated when a craving occurs, would be benefi cial to someone attempting to abstain from smoking. Furthermore, a system that can be tailored to an individual's smoking schedule would provide a more effi cacious treatment than a standard "one size fi ts all" therapy.
Similarly, a multiple rate delivery system would be benefi cial to patients and medical staff when treating opiate withdrawal symptoms using clonidine. Currently, clonidine is used in rehabilitation clinics to mitigate the symptoms of opiate withdrawal including rhinorrhea, lacrimation, perspiration, piloerection, gastrointestinal complaints, nausea, yawning, sneezing, muscle pain, tachycardia and hypertension. 4 Clonidine is a partial pre-and postsynaptic alpha-adrenoreceptor agonist that produces a reduction in peripheral sympathetic outflow. Patients experience acute withdrawal symptoms between 6 to 12 hours after opioid abstinence with peak effects being felt between 36 to 72 hours and gradual abatement over 5 to 10 days. 3,9 For this reason patients are typically given multiple oral clonidine doses, or in a few cases a combination of transdermal and oral drug, per day during treatment. 4,5 A typical fi ve day oral treatment profi le can be seen in Figure 1b. Clearly, repetitive dosing with variable amounts each day puts a strain on the patient and clinical staff. Patients who are experiencing withdrawal are likely irritable and unreceptive to medical staff who must administer the multiple doses of a drug per day. Additionally, patients who suffer nausea and gastrointestinal problems during withdrawal experience diffi culties in maintaining proper doses due to vomiting. A transdermal system that can manipulate its delivery rate to provide the required dose on specifi c days during a typical clinic stay would not only eliminate the physiological dosing issues, but also ease the burden of administering, monitoring, and tailoring a drug regimen to each individual patient.
Transdermal patches tend to be one of the most desired therapies by patients and health care providers, not only for easy and pain-free administration, but also for several pharmacokinetic and pharmacodynamic benefi ts. 10,11 Traditional transdermal delivery allows one time application to constantly administer drug over a given time, thereby eliminating the hassle and diffi culties associated with oral or intravenous maintenance multiple dosing. This "worry free" administration increases patient compliance, which is essential when trying to combat drug addiction. Additionally, once a drug passes through the upper echelons of the skin, it is directly taken into systemic circulation by the highly vascularized dermal layers. This affords the opportunity to lower the amount of drug required at the administration site, because of the avoidance of fi rst pass metabolism. Not only is this more cost effective, but also safer for the patient as high levels of compounds can increase the chance for toxicity. Furthermore, as mentioned previously, gastrointestinal complications can lead to diffi cult and improper oral dosing when treating patients experiencing opioid withdrawal symptoms. Such symptoms would not hamper therapeutic regimen designs when delivering the desired pharmaceutical agent transdermally. Despite the attractiveness of the transdermal route, one major shortcoming exists in effi caciously treating smoking and opioid addiction with currently available patches; none of the current marketable transdermal patches can deliver the required variable dosing rate of drug to effectively treat these addictions in a single, self-contained system. The use of CNT membranes in a transdermal system offers a feasible solution to this variable rate drug delivery predicament.

Carbon Nanotube Membranes
The demand for development of highly functionalized chemical screening, separating, purifying and delivery devices has risen in the past several decades as vast advances in the areas of biology, chemistry, electronics and medicine continue to be made. 12 Nanoporous membranes consisting of aligned, multi-walled CNTs in a polymer substrate, possessing distinctive properties have demonstrated their applicability in the areas of electronic devices, chemical separation, biosensors, salt water desalination and purifi cation, and most recently, in the area of drug delivery. More specifi cally, three signifi cant properties unique to CNT membranes that can be utilized to conquer the diffi culties of a variable rate drug therapy are 1) extremely fast fl uid fl ow due to graphite planes and large van der Waals distance, 13 2) well defi ned geometry of

Day
Oral Dose ligands at pore entrance resulting from membrane synthesis, 14 and 3) the concentration of an electric fi eld at the tips of CNTs due to their electrically conductive properties while surrounded by an insulating polymer support. 15 Synthesis of these membranes has been described previously, 16 but briefly carbon nanotube arrays are grown by chemical vapor deposition (CVD) on a quartz substrate using specifi c metal catalysts to control the diameter of nanotubes. The spaces between the individual tubes are then embedded with a polymer substrate and the quartz support etched off using hydrofl uoric acid. H 2 O plasma oxidation is utilized to remove excess polymer on the surface of the fi lm as well as the initial growth catalyst particle. This process results in a robust membrane with open-ended CNTs possessing carboxylic functionality at the tips. These tips can be further functionalized with various molecules using carbodiimide chemistry to create charged molecular chains that, with proper electrical infl uence, will work as a gate to chemical transport through the nanotubes. 15 More specifi cally and as seen in Figure 2, by applying a small electrical bias (Ͻ50 mV) to the membrane, the end groups, which are negatively charged, are either attracted to the nanotubes themselves causing blockage of the pore and thus no transport, or repelled from one another and surrounding nanotubes causing the functional group to stand on end allowing fl uid passage and drug delivery through the pore. By incorporating this voltage gated membrane into a transdermal system (Fig. 3a), programmable controlled dosing either by physician or patient is attainable.
Recall that the major set back in current NRT therapies is the inability to mimic the actual pulsatile nicotine delivery profi le and resultant peak plasma levels of a smoker. With the ability of a patch to deliver a high dose of nicotine and then safely shut off to return to baseline plasma levels throughout the day, a smoker would be able to satisfy their nicotine craving and thereby increase their chances of becoming smoking independent. Furthermore, these patches can be patient controlled through a push activated button pressed at the initiation of a patient's craving, which would apply the electrical bias to the membrane allowing drug fl ow to the skin, eliminating patient to patient tailoring issues seen with conventional therapies. As with all drug addiction and subsequent pharmacological replacement treatment, the patient must be motivated to quit to truly gain independence from whatever drug one is dependent upon.
Similarly, the multiple and variable dosages that are utilized to properly treat opioid withdrawal symptoms can be simplifi ed to a one time application and computer chip controlled transdermal fl ux change system. Using multiple membranes within the patch will allow for a therapy that is patient tailorable and highly desirable for both patient and medical staff. Using the setup seen in 3b, the typical clonidine dosing profi le previously discussed can be achieved by activating one side of the patch the fi rst and fourth day, activating both membranes the second and third day for the highest transdermal fl ux, and by activating only the second membrane the fi fth day to taper down the dose.
As alluded to earlier, abusers going through withdrawal may experience increased sweating, variable circulation and based on the degree of withdrawal, craving, or intoxication may be quite different from patient to patient. 17 This variability in skin conductance could introduce problems in traditional transdermal delivery; however one of the main advantages of the CNT membrane system is its increased dosing fl exibility. More specifi cally, the system cannot only be tailored to an individual patient's needs, but also possesses the added benefi t of fl exibility in drug delivery as conditions change over the course of treatment. Because of the high doses of nicotine required to effectively treat smoking addiction with this nanotube system and the fact that the system will be patient activated, the potential for lethal overdose is a possibility without a proper design. For this reason the patch will require a multi-use lockout system; the fi rst facet to automatically shut off the patch when the proper nicotine amount is delivered, the second to prevent the addicted individual from activating the patch too often and delivering too much nicotine in a short time period. The proper amount of nicotine refers to the required amount of drug to be delivered transdermally that corresponds to the nicotine plasma levels seen immediately following the smoking of a cigarette.Additionally, as discussed, this patch setup will involve the use of a small electric current to deliver the drug. Currently marketed transdermal iontophoresis systems also use electric current to deliver drug. However, unlike iontophoresis, the electrodes used to deliver the bias are never in contact with the skin, nor is the voltage required as high; 50 mV for the CNT system compared to 10 V for iontophoresis. 15,18 For these reasons the skin irritation and lesions that are often seen with the iontophoretic set up are not expected when using CNT membranes in a transdermal patch. 19 The use of CNT membranes in transdermal therapy is fairly cost effective as the membranes themselves can be fashioned for as little as $0.60/m 2 . The bulk of the cost of this transdermal set up will come from cost of drug and the bias control system. The cost of drug will obviously be comparable to that of current nicotine and clonidine transdermal therapies. The voltage required to activate and deactivate the system is no more than that which can be supplied by a watch battery, and the mechanics of the system, electrodes to provide the bias, is comparable to that of the new patient activated iontophoteric fentanyl delivery system marketed as IONSYS TM . 20 Additionally, the voltage system would be able to be reused, making the drug the only continuous expense throughout the therapy.

Therapeutic Feasibility
Despite the simplicity of the system, it was still essential to demonstrate that the skin alone could permit the required maximum transdermal fl ux, proper rate changes, and especially in the case of nicotine delivery, do so with minimal lag time. A series of proof-of-concept in vitro diffusion studies were conducted on excised human abdominal skin in order to determine if the three parameters discussed were feasible for each drug. The average smoker has a steady-state nicotine plasma concentration of about 10-20 ng/ml with peak levels reaching approximately 54 ng/ml throughout the day with each cigarette. 7 Using equation 1 in which J is the steady-state transdermal fl ux, C ss is the steady-state plasma concentration, Cl is the total body clearance of nicotine, and A is the diffusional area (7 cm 2 ), these plasma levels translate to a required transdermal fl ux of 0.8 μmol/cm 2 /h and 4 μmol/cm 2 /h for steady-state and peak nicotine concentrations, respectively. For a traditional fi ve-day opioid withdrawal symptom treatment regimen using clonidine, the fi rst and fourth days require a total of 0.2 mg per day of oral clonidine, days two and three a total of 0.3 mg per day, and the fi fth day 0.1 mg of oral clonidine as seen in Figure 1b. 5 Using a patch size of 10 cm 2 , 100% bioavailability, 21 and each day's oral dose over a 24 hour period, effective transdermal fl uxes of 3.5 nmol/cm 2 /h, 5.4 nmol/cm 2 /h, and 1.7 nmol/cm 2 /h can be estimated.
Human skin from abdominoplasty surgery was obtained from the National Cancer Institute's Cooperative Human Tissue Network (CHTN). The samples were dermatomed immediately upon arrival to a thickness of approximately 200 μm and frozen at −20 °C. On the day of the experiment, the skin was thawed and immediately used for the diffusion studies. The abdominal skin was loaded into a cell of a PermeGear ® In-Line Diffusion apparatus (Bethlehem, PA, U.S.A) and charged with 250 μl of desired drug solution, and temperature was maintained at 35 °C (skin surface temperature) by a circulating water bath. 22 Solutions of 0.9% saline for nicotine and 0.9% saline with 10% ethanol for clonidine were prepared for use as the receiver solutions for diffusion. Multiple donor concentrations were used in order to simulate a CNT change in drug supply to the skin surface to determine whether a fl ux change would be observed through the skin, and to see its associated transition (lag) time for rate change. For example, after it was determined that a nicotine donor concentration of 500 mg/ml would produce the desired 4 μmol/cm 2 /h fl ux, this donor solution was applied to the skin and the receiver solution was sampled over 30 minutes. Next, the donor solution was removed and the skin and diffusion cell donor chamber quickly rinsed with water to remove any residual nicotine on the surface, and a much lower donor solution concentration of 50 mg/ml was used in order to observe a fl ux change and lag time to fl ux change. After another 30 minutes of receiver solution sampling, the skin and cell chamber were washed again and a donor solution of 500 mg/ml was reapplied to the skin, and the receiver solution sampled again. This process mimics the activation of the carbon nanotube membrane to allow more drug fl ow, followed by deactivation and then reactivation, much like a person smoking two cigarettes over an hour and a half period. Transition times were observed and results can be seen in Table 1. The term transition time for this work refers to the time it takes between a change in drug rate supply to reach the new desired steady-state concentration in the receiver solution.
It was also important to investigate the rate of drug delivery out of the skin after a donor concentration was removed, in other words the depletion of nicotine from the skin after the patch is turned off. The residual drug amount that is in the skin after the patch is shut off or removed will enter systemic circulation and must be taken into consideration when determining the proper time to allot between patient initiated activation. For this study, varying donor solution concentrations (200, 300, 400, 500 mg/ml) were applied to the skin long enough to reach steady state, and then removed and the cell chamber and skin briefl y washed with water. Sampling of the receiver solution was continued for approximately two hours with no donor solution on the skin in order to observe nicotine depletion. These studies resulted Table 1. Therapeutic transdermal fl uxes and transition times for nicotine and clonidine (n = 4).

Concentration
(mg/ml) in the model seen in Equation 2, where the variables correspond to Figure 4 and k is the time constant which was 0.088 ± 0.0092 min −1 (Graphpad Prism®, La Jolla, CA, U.S.A) and independent of donor concentration (p Ͼ 0.05) (one-way ANOVA, SigmaStat3.5 © , San Jose, CA, U.S.A). This indicated that the skin will not have a signifi cant depot effect once the patch is removed or the system is completely shut off and that proper lockout periods can be programmed based on this equation.
Clonidine also had to be investigated to ensure that the proper therapeutic levels of drug could be reached when delivered transdermally, in addition to if the skin responds to a change in donor concentration allowing a subsequent change in transdermal fl ux. A similar process of applying a donor solution followed by a brief washing of the skin and diffusion cell donor compartment was repeated, except that this time each donor solution remained on the skin for 24 hours thereby more closely representing therapeutic clonidine delivery. The required therapeutic fl ux values for each of the different doses were achieved and can also be seen in Table 1, but it should be noted a short transition time between doses is not as important as it is with nicotine, as tapering doses over the 5-10 day treatment are currently used in the clinical setting.

Future CNT Membrane Application
It is easy to envision the use of the CNT membrane in therapies not discussed here. Its applicability is wide and the ability to change the rate of delivery will prove highly benefi cial in many therapeutic areas. Any skin permeable drug is suitable for CNT membrane delivery as the membrane is not rate limiting in the way the skin can be based on size of the molecule or lipophilicity of the drug. CNT membranes do not inhibit the use of transdermal techniques to increase the breadth of drugs available for use in transdermal systems such as prodrugs and codrugs, or the use of microneedles prior to patch application. [23][24][25][26] Using these techniques in conjunction with CNT membranes would not only allow a wider range of drugs to be delivered, but also provide more fl exibility in programmable rates. By increasing the permeation through the skin, the range of high to low fl ux through the membrane is increased allowing fi ner patient to patient tailoring. Additionally, a multiple drug delivery system could be developed in which two drug reservoirs of different content and two separate CNT membranes are utilized. This would allow for simultaneous delivery of two drugs at varying rates, or for transfer between medications without removal and application of a new transdermal system. Obviously, size and cost of the device would increase as the system becomes more advanced, however the benefi ts of variable delivery would outweigh these disadvantages in some therapeutic applications.

Conclusions
Carbon nanotube membranes present a feasible solution to pitfalls in current smoking cessation and opioid withdrawal symptom treatments. Therapeutic transdermal fl uxes can be achieved and the skin responds well to changes in donor concentrations indicating that the skin will not be a limiting factor when bias is applied to change the fl ux from the patch. This novel method to create a rate changing drug delivery device will not only increase patient compliance because of the multiple benefi ts of using the transdermal route, but also provide more effective therapy and improve the quality of life for many addicted individuals.

Acknowledgment
This work was supported by NIH grant R01DA18822.

Disclosure
The authors report no confl icts of interest.


extfieldsofstudy: Biology Medicine

text:
Remarkable sequence similarity between the dinoflagellate-infecting marine girus and the terrestrial pathogen African swine fever virus

Heterocapsa circularisquama DNA virus (HcDNAV; previously designated as HcV) is a giant virus (girus) with a ~356-kbp double-stranded DNA (dsDNA) genome. HcDNAV lytically infects the bivalve-killing marine dinoflagellate H. circularisquama, and currently represents the sole DNA virus isolated from dinoflagellates, one of the most abundant protists in marine ecosystems. Its morphological features, genome type, and host range previously suggested that HcDNAV might be a member of the family Phycodnaviridae of Nucleo-Cytoplasmic Large DNA Viruses (NCLDVs), though no supporting sequence data was available. NCLDVs currently include two families found in aquatic environments (Phycodnaviridae, Mimiviridae), one mostly infecting terrestrial animals (Poxviridae), another isolated from fish, amphibians and insects (Iridoviridae), and the last one (Asfarviridae) exclusively represented by the animal pathogen African swine fever virus (ASFV), the agent of a fatal hemorrhagic disease in domestic swine. In this study, we determined the complete sequence of the type B DNA polymerase (PolB) gene of HcDNAV. The viral PolB was transcribed at least from 6 h post inoculation (hpi), suggesting its crucial function for viral replication. Most unexpectedly, the HcDNAV PolB sequence was found to be closely related to the PolB sequence of ASFV. In addition, the amino acid sequence of HcDNAV PolB showed a rare amino acid substitution within a motif containing highly conserved motif: YSDTDS was found in HcDNAV PolB instead of YGDTDS in most dsDNA viruses. Together with the previous observation of ASFV-like sequences in the Sorcerer II Global Ocean Sampling metagenomic datasets, our results further reinforce the ideas that the terrestrial ASFV has its evolutionary origin in marine environments.

synthesis, some are heterotrophic using endocytotic feeding, and many dinoflagellates are mixotrophic having both modes of nutrition. Blooms of certain photosynthetic dinoflagellates kill fish and bivalves, or pollute shellfishes for food with particular toxins, and can lead to serious economic damages in aquaculture [3,4]. Heterocapsa circularisquama forms blooms causing massive death of shellfish such as pearl oysters and mussels, and is one of the most intensively studied dinoflagellate species [5].
HcDNAV is a marine giant virus (or "girus" [6,7]) containing dsDNA genome, and lytically infects H. circularisquama [8,9]. HcDNAV is considered to play a significant role in the demise of H. circularisquama blooms [9,10]. HcDNAV has a large icosahedral capsid (180-210 nm in diameter), which packs a ~356-kbp genome [8,11]. During its multiplication, virions emerge from a specific cytoplasm compartment, called "viroplasm", which is created by the virus [9]. HcDNAV is the sole DNA virus currently isolated from dinoflagellates, and to our knowledge, is the only DNA virus isolated from the superphylum Alveolata [12]. Based on its host range, genome type/size and microscopic features, HcDNAV was previously suggested to be a member of Phycodnaviridae [13]. However, there has been no molecular data supporting this tentative classification.
Phycodnaviridae includes intensively-studied algal virus members such as chlorella viruses and Emiliania huxleyi viruses [14][15][16][17], and belongs to a larger group of eukaryotic DNA viruses called NCLDVs [18]. NCLDVs complete their replication cycle within the host cytoplasm, and share an array of conserved core genes for transcription, RNA processing, replication, DNA packaging, and structural components. Other viral families of NCLDVs are Mimiviridae, Poxviridae, Iridoviridae, and Asfarviridae. Mimiviridae is represented by the freshwater amoebainfecting mimivirus [19] and its close relative mamavirus [20]. Based on the sequences of PolB, the most conserved NCLDV core genes, three algal viruses have been suggested to belong to Mimiviridae [21]. Poxviridae include a number of successful pathogens known to infect a tremendous variety of terrestrial animals, such as insects, reptiles, birds, and mammals [22]. Iridoviruses infect invertebrate and cold-blooded vertebrate hosts, and includes numerous emerging pathogens of fishes and amphibians [23]. The last family Asfarviridae [24,25] is currently represented by a sole species, African swine fever virus (ASFV) with a 170 kbp dsDNA genome [26]. ASFV is a large (~200 nm in diameter), intracytoplasmically-replicating arbovirus, naturally maintained in a sylvatic cycle between wild swine (warthogs and bushpigs) and argasid ticks (Ornithodoros). In these hosts, ASFV infection is usually asymptomatic [27]. However, ASFV causes an acute hemorrhagic infection in domestic swine with mortality rates up to 100% for some viral isolates.
In an attempt to further characterize HcDNAV, we performed a low coverage shotgun sequencing of its genome. Specifically, from 4 liters of HcDNAV suspension (lysate of HcDNAVinfected H. circularisquama on 6 dpi), virus particles were collected as described in [11]. The viral genomic DNA was purified in a PFGE-gel and was subjected to shotgun sequencing (coverage = 0.11 X). Resulting sequence reads covered part of the region containing a PolB-like sequence. With the use of tail-PCR method [28], we successfully determined a 5,800 bp sequence (DDBJ accession number AB522601) containing an open reading frame (ORF) for the complete HcDNAV PolB gene. By means of a reverse transcription-PCR (RT-PCR) experiment, the PolB gene was shown to be transcribed to mRNA (additional file 1); thus, it is most likely crucial for the replication of HcDNAV.
HcDNAV PolB gene was found to be 3,675 bp long (forward strand, position = nt 1,913-5,590 in AB522601), punctuated by normal start and stop codons, and no intron or intein-like sequence was observed. The predicted protein product is 1,225 amino acids (aa) long. Unexpectedly, the translated amino acid sequence showed the closest BLASTP hits against PolB sequences from different ASFV isolates, with the best homolog being DPOL_ASFL6 (identity = 27%, bit score = 311, E-value = 4.10E-82) in the NCBI non-redundant sequence database. The best non-ASFV hit corresponded to the PolB sequence of Pyramimonas orientalis virus (DPOL_POV01, identity = 23%, bit score = 131, E-value = 4.10E-28). A multiple sequence alignment of the HcDNAV PolB and its close homologs confirmed the presence of conserved residues for exonuclease and polymerase activities [29] (additional file 2). Curiously, the HcDNAV PolB sequence exhibited a rarely observed amino acid substitution within the motif containing two highly conserved metal binding aspartic acid residues; HcDNAV exhibits the motif YSDTDSinstead of the YGDTDS-sequence usually found in dsDNA viruses. In addition, we identified two ORFs in the upstream region of the PolB ORF in a divergent orientation. Their products were respectively predicted to be 245 and 194 aa in length (positions = nt 463-1,200 and 1,255-1,839). The former showed a significant similarity to HNH endonucleases with its BLASTP best hit to mimivirus L245 (YP_142599, E-value = 4E-11); the latter showed a significant similarity to hypothetical proteins from NCLDVs with its best hit to mimivirus R325 (annotated as a metal-dependent hydrolase, YP_142679.1, E-value = 1E-12). Incidentally, R325 is located near the PolB gene (R322) in the mimivirus genome [30].
To examine the unexpected sequence similarity between the HcDNAV and ASFV PolBs, we conducted a series of maximum likelihood phylogenetic analyses. First, we aligned the HcDNAV PolB sequence with its homologs from NCLDVs. A phylogenetic tree based on the 362 amino acid residue sites from the alignment supported the monophyletic grouping of HcDNAV and ASFV with a 100% bootstrap value (Fig. 1). The grouping of each of the Maximum likelihood tree of PolB amino acid sequences from NCLDVs Figure 1 Maximum likelihood tree of PolB amino acid sequences from NCLDVs. Alignment was constructed with the use of T-Coffee. All the gap-containing amino acid residue sites were removed before tree construction. The phylogenetic tree was constructed using PhyML [38] available at Phylogeny.fr [39] using WAG matrix and gamma distribution. Branch labels indicate bootstrap percentages (≥ 50%) after 100 replicates. The tree is essentially an unrooted tree, albeit mid-point rooted only for presentation purpose. The same method was used for the phylogenetic trees in Fig. 2, Fig. 3  other four NCLDV families was also supported by a high bootstrap value (100% for Iridoviridae, 81% for Phycodnaviridae, 90% for Mimiviridae and 100% for Poxviridae). Next, we used a wider range of viral homologs including those of bacteriophages. The resulting tree based on 320 amino acid residues again supported the grouping of HcDNAV/ASFV with a 98% bootstrap value (Fig. 2).
In addition, we obtained a short sequence partially corresponding to an RNA polymerase II large subunit gene from HcDNAV genomic DNA (AB522602), for which we obtained a similar result. The 892 bp sequence showed BLASTX best hit against ASFV RNA polymerase sequence (RPB1_ASFM2, E-value = 2E-12). A monophyletic grouping between the HcDNAV sequence (97 aa) and the ASFV

Figure 3 Maximum likelihood tree of PolB amino acid sequences from NCLDVs and several sequences from environmental samples (indicated by open diamond marks).
HcDNAV and ASFV sequences are indicated by filled diamond marks.
RNA polymerase sequence was again received a high bootstrap value of 87% (additional file 3).
Our homology search and phylogenetic analyses thus confirm that the newly determined HcDNAV sequences are most closely related to their ASFV homologs. This result is in clear contradiction with the previous proposal that HcDNAV may belong to the Phycodnaviridae [13].
A previous "phylogenetic mapping" survey of the metagenomic sequence data sets generated by the Global Ocean Sampling (GOS) expedition [31] revealed several PolBlike sequences most closely related to the PolB sequence of ASFV [32]. This observation suggested the presence of ASFV-related viruses in marine environments. In order to examine whether the "ASFV-like" marine PolB sequences were close to the HcDNAV PolB sequence, we retrieved 267 sequences from the environmental sequence collection of NCBI/GenBank using the PolB sequences of HcD-NAV and ASFV as queries (E-value < 1E-10). These environmental sequences were in turn searched against the NCBI non-redundant sequence database and the HcDNAV PolB sequence. Of the 267 sequences, 15 showed their best hit to the ASFV PolB, one showed its best hit to HcDNAV (gi|136563424), and the remaining sequences had their best hit to other viruses or cellular organisms. Therefore, most of the ASFV/HcDNAV-like PolB sequences in the marine environmental collection are more closely related to the ASFV PolB than to the HcD-NAV homolog. A phylogenetic tree using several environmental sequences supported their grouping with the terrestrial ASFV PolB (bootstrap value = 84%, Fig. 3).
PolB is one of the most reliable phylogenetic markers for large eukaryotic DNA viruses [32,33]. The fact that the HcDNAV PolB was not grouped with the PolBs from phycodnaviruses strongly argues against the previous tentative classification of HcDNAV in the Phycodnaviridae family [13]. It is clear that the definitive classification of HcDNAV will require the complete sequencing of its genome. It may also turn out that the HcDNAV genome corresponds to a mosaic of NCLDV genes with different evolutionary histories, precluding a simple classification scheme. Pending its complete genome sequencing, we recently proposed to the ICTV to create a new genus "Dinodnavirus" where to tentatively classify the HcDNAV.
Our finding now establishes an evolutionary link between a terrestrial pathogen and a marine girus. A recent metagenomic analysis of corals provided evidence for the existence of viruses related to herpesviruses [34], which have been mostly isolated as pathogens of terrestrial animals. So far, giruses of 7 algal classes [12,35] have been isolated; still, we know next to nothing about viruses infecting other protists in aquatic environments. Given the huge diversity of protists [36,37], a comparable diversity probably exists for marine viruses living in these environments. Exploring this hidden viral world is necessary to our understanding of the evolutionary relationships between aquatic viruses and their terrestrial relatives.

Competing interests
The authors declare that they have no competing interests.

Authors' contributions
NK conceived the study. YS, KT and NN conducted purification and sequencing of HcDNAV and RT-PCR experiment. HO designed and carried out bioinformatics analyses. HO, JMC, YT and NK contributed to the interpretation of data and wrote the manuscript. All authors read and approved the final manuscript.

Additional file 2
Conserve blocks from the multiple sequence alignment of PolB sequences from NCLDVs. The data provided shows the presence of conserved residues for exonuclease and polymerase activities in the HcDNAV PolB and its close homologs. Species abbreviation is followed by a database sequence identifier. Intein sequences were removed from the sequences prior to alignment. The alignment was generated by T-Coffee [40] and ClustalX [41]. AmEPV