BACKGROUND
interactions between eu- and prokaryotic cells are frequent, multifaceted events ranging from symbiotic synergy such as symbiotic nitrogen fixation in legumes or fermentation by gastrointestinal bacteria to pathogenic interference, for instance, in the course of salmonellosis.
this interplay of organisms requires mutual signaling mechanisms and a continuous adaptation of the metabolism of the involved cells to varying environmental conditions.
consequently, programmed expression patterns have to be induced to continuously readjust the proteome and metabolome of both cell types.
the characterization of corresponding time-dependent expression patterns allows for a deeper understanding of the underlying molecular processes, and was the focus of numerous studies, but until recently gene expression profiling emphasized either the host cell or the prokaryotic transcriptome  <cit> .
salmonella represents a genus of gram-negative and facultative anaerobic enterobacteria and is closely related to the genus escherichia.
the human-virulent pathogen is a model organism for characterization of host-pathogen interactions, and is associated with a variety of diseases including gastroenteritis and enteric fever.
one of the first attempts to profile genome-wide expression changes in host cells during interaction of pro- and eukaryotes was carried out with s. enterica and human intestinal epithelial cells  <cit> .
conversely, the transcriptome of s. enterica became subject to several studies of host-pathogen interactions after completion of the genome sequences of s. enterica serotype typhi ct <dig>  <cit>  and serotype typhimurium lt <dig>  <cit> , which complemented some of the previously determined host responses  <cit> .
in the meantime, next-generation sequencing -coupled transcription profiling techniques emerged as the principal tools to interrogate gene expression, and especially whole transcriptome shotgun sequencing  has considerably contributed to our understanding of prokaryotic transcriptomes  <cit> .
nonetheless, simultaneous transcription profiling without prior disruption of the interaction remains technically challenging, and thus characterization of disease-related expression patterns in interacting eu- and prokaryotic cells is inevitably linked to comprehensive sequencing efforts  <cit> .
here we present dual 3’seq, a tag-based, ngs-coupled method that allows for simultaneous transcription profiling of interacting pro- and eukaryotes without physical separation of the interacting cells.
compared to rna-seq, the reduction in complexity of tag-based approaches significantly decreases the required sequencing depth for a good coverage of both the pro- and eukaryotic transcriptomes  <cit> , which is a prerequisite for profiling of low abundant pathogen-derived transcripts.
additionally, only a single tag is generated out of each transcript, which facilitates unequivocal quantification of reads from a specific rna without sacrificing qualitative information of pathogen-derived transcripts, since prokaryotes lack alternative splicing events  <cit> .
deepsupersage  and mace  represent two established ngs-coupled transcriptome profiling techniques that generate exactly one tag out of the 3′ end of every transcript.
while deepsupersage yields a  <dig> nucleotide tag that is specifically located within the 3′ end depending on the presence of an according restriction site for the employed anchoring enzyme, mace generates randomly distributed tags out of the last hundreds of bases.
mace consequently allows for preparation of libraries with varying read lengths and provides additional information regarding transcripts that do not possess an according restriction site for anchoring.
in order to assess the respective efficiencies in transcriptome profiling of cultivated and interacting s. enterica typhimurium and human host cells, we combined the dual 3’seq approach with both protocols.
human epithelial cells  were infected with the invasive pathogen s. enterica typhimurium sl <dig> , and interacting cells were screened for differentially expressed transcripts at several points of time post infection to provide an overview of the transcriptional processes during invasion of the intestinal epithelium as one of the first steps in emerging salmonellosis.
the combination of the published sl <dig> transcriptome  <cit>  with the operon structure identified by differential rna-seq  allowed for accurate quantification of polycistronically transcribed genes from the prokaryote, and the corresponding expression profiles provide a basis for time-dependent analysis of disease-related transcripts, despite the fact that transcripts from the pathogen are of extremely low abundance during interaction.
RESULTS
dual 3’seq of interacting pro- and eukaryotes
in order to reduce the sequencing depth for simultaneous characterization of the transcriptomes of interacting pro- and eukaryotic cells without prior separation, we established the new transcriptome profiling technique dual 3’seq .
we employed this technique to identify time-dependent transcription patterns of interacting hela-s <dig> and sl <dig> cells subsequent to infection of the host cells.
total rna from cultured and interacting cells was isolated after different points of time to interrogate the molecular events during early, mid-level and late interaction stages of infection .
the first interaction point was chosen as early as technically possible to capture the initial pathogen responses after invasion of the host cells, while mid-level and late interaction allow for investigation of the early and late host cell response as well as intracellular pathogen replication.figure  <dig> scheme of dual 3’seq library preparation and bioinformatic processing of the generated sequencing data.
total rna was size-selected  subsequent to dnase i digestion of remaining dna in the isolate.
following rrna depletion , the rna was split into the poly+ and poly− fraction by oligo capture to separate the polyadenylated and functional mrnas of eukaryotic cells from the non-polyadenylated transcripts that represent the functional transcriptome of prokaryotes.
ensuing in-vitro polyadenylation of the poly− fraction, both fractions were subjected to oligo-based reverse transcription.
the generated cdna was fragmented according to two established 3′ transcriptome profiling techniques.
deepsupersage tags were generated via cleavage of rnas by the anchoring enzyme nlaiii and subsequent digestion using ecop15i, while mace involved random fragmentation for generation of tags.
3′ fragments were enriched by binding to a streptavidin matrix and ligated to a sequencing adaptor.
adaptor-ligated fragments were pcr-amplified using genxpro’s truequant technology for pcr-bias free amplification, page-purified, and sequenced on the illumina hiseq <dig> platform.
barcoded reads were allocated to their respective library, filtered for pcr-derived reads, and trimmed for high-quality sequences.
afterwards, reads were annotated to a combined reference comprising the transcriptome and genome sequences of sl <dig> and human host cells in a multi-step procedure.
reads uniquely mapped to one of both organisms were combined to three distinct expression matrices for functional analysis of the poly− transcriptome from pathogen and host cells as well as the poly+ fraction of the host cells.
for each expression matrix, annotated reads were quantified and median-normalized using deseq, followed by pair-wise, time-dependent comparison of the different interaction stages.
statistical significance was subsequently corrected for multiple testing according to benjamini and hochberg.
the generated sequencing data from the poly− fractions comprise between  <dig>  and  <dig>  million high-quality and pcr duplicate-free reads, while the number of reads in the corresponding poly+ libraries ranges from almost  <dig> million down to ~ <dig>  reads in the sl <dig> library .
the relatively low number of filtered poly+ reads from cultivated and especially interacting sl <dig> cells reflects post-transcriptional regulation in prokaryotes, where polyadenylation facilitates degradation of rna  <cit> .
additionally, the low abundance of polyadenylated intermediates from the prokaryote is corroborated by the low ratios of sl1344-derived transcripts in all interaction stages , which impedes a functional analysis of this fraction.
prokaryotic reads in libraries of the poly− fraction are about ten times more abundant .
compared to the ratios of the corresponding poly+ libraries, non-polyadenylated sl <dig> reads are about 160-fold, 50-fold and 10-fold higher from early to late interaction, respectively.
this increase in polyadenylated reads during interaction suggests an increased turnover of transcripts probably caused by the growing number of  pathogen cells in the confined environment of hela cell culture.
discarding ambiguously annotated poly− reads that aligned to more than one gene,  <dig>  or 89% of the  <dig>  known loci in sl <dig> are transcribed in cultivated or interacting cells, while  <dig>  or 75% of  <dig>  protein-coding loci are represented in the polyadenylated transcriptomes of human host cells.
over 95% of the transcripts in human poly+ libraries encode proteins, and mrnas encoding ribosomal proteins account for a large proportion of these .
with an mrna content ranging from 13% to 20% protein-coding transcripts are much less abundant in the corresponding poly− libraries, which mostly comprise small nucleolar rnas , small nuclear rnas  and other non-coding rnas .figure  <dig> annotation statistics for the generated expression matrices and differentially expressed genes of interacting pathogen and host cells.
read numbers annotated to human  and pathogen  cells or both  for all sequenced libraries, respectively.
the bar graph shows log10-transformed numbers, while the original number of reads is depicted in a circle graph.
classification of unambiguously annotated transcripts from host cells into mrnas , non-coding rnas , bivalent transcripts encoding mrnas and non-coding rnas , long non-coding rnas , histone-encoding mrnas , and other transcript classes  for the poly+ and poly− transcriptome, respectively.
mrnas encoding cytoplasmic  or mitochondrial  ribosomal proteins are depicted individually for the poly+ and combined  for the poly− transcriptome that additionally comprises small nucleolar rnas  and small nuclear rnas .
respective abundances are shown for non-interacting  as well as interacting  cells.
the different interaction stages  prepared with 3’seq coupled to deepsupersage  are depicted from inside outwards along with the early interaction library that was prepared using mace .
venn diagram showing the number of commonly expressed genes for all combinations of libraries corresponding to one of the expression matrices.
ambiguously annotated reads are not included.
pair-wise comparison of differentially expressed transcripts within the expression matrices  along with the number of exclusively expressed transcripts  for each library.
listed are unambiguously annotated reads with a log <dig> fold change stronger than | <dig> | and an fdr-corrected p-value below  <dig> .
up- and downregulated transcripts are indicated by corresponding arrows for the library on the ordinate, respectively.
differential gene expression during interaction
the intersections of unambiguously annotated rnas from the poly− transcriptome of sl <dig> cells as well as both transcriptome fractions of human host cells reveal distinct overlaps of the expressed genes .
the four sl <dig> libraries comprise ~ <dig> commonly expressed genes, while non-interacting and late interacting cells share the highest number with more than  <dig> .
unsurprisingly, cells from early and mid-level interaction express fewer genes in common with non-interacting cells than those from the late interaction stage when infection is well established.
conversely, the number of commonly expressed genes in the poly+ transcriptome of host cells is highest between early and mid-level interaction .
with  <dig>  ±  <dig>  genes, the numbers of overlapping genes in the poly− fraction of host cells spread notably less than in the poly+ fraction, which suggests a more diverse effect of infection on the polyadenylated transcriptome of host cells.
although the poly− fraction of host cells generally comprises more significantly up and downregulated transcripts than the poly+ fraction , expression of polyadenylated transcripts during infection is more affected, since mrna degradation naturally leads to corresponding changes in the non-polyadenylated transcriptome.
in fact, more than 80% of the non-polyadenylated, significantly differentially expressed transcripts in all the pairwise comparisons represent protein-coding mrnas, even though these rnas only constitute up to a fifth of their respective poly− fractions in total .
expression of stm <dig> and stm <dig> is non-detectable except for the transcriptome of early interacting sl <dig> cells, while the transcripts encoded by gst, ydin and stm <dig> are only present in mid-level interaction .
stm <dig> is part of salmonella pathogenicity island  <dig> , and substantially contributes to bacterial fitness by promoting transcription of spi- <dig> together with ssrb  <cit> .
stm <dig> and stm <dig> are both encoded within gifsy- <dig>  a lambdoid prophage that promotes virulence in addition to the remnant phage encoded by spi- <dig>  <cit> .
the time-dependent differential expression of all sequenced genes from the prokaryote reveals characteristic patterns for each interaction stage .
the expression of many spi-encoded and phage-derived transcripts from the chromosome as well as plasmid-derived transcripts is induced or repressed in the course of infection, and analysis of these disease-related expression patterns allows for dissolving the time-dependent components of sl <dig> gene expression during different stages of infection.figure  <dig> differential expression of interacting pathogen cells.
circos plot of sl <dig> gene expression for chromosomally-encoded genes as well as plasmid-derived transcripts.
the respective coordinates are given in megabases  along with the corresponding log <dig> fold changes of non-interacting and early interacting, early and mid-level interacting, as well as mid-level and late interacting cells .
salmonella pathogenicity islands  are highlighted in orange, while other important loci are shown in green.
detailed time-dependent expression profiles of selected sl <dig> genes.
the median-normalized transcript abundance is shown for non-interacting cells along with those from the different interaction stages.
ptsn is transcribed as polycistronic mrna from the chromosome that comprises ten other genes .
the mrna encoding pefd and the polycistronic mrna coding for repa, tap, and repa <dig> are both transcribed from intracellular plasmid psltsl <dig>   k-means clustering of the time-dependent fold changes from plasmid-encoded transcripts.
the corresponding log <dig> fold changes were assigned to one of four clusters based on their time-dependent expression patterns.
the closest centroid  is shown along with the assigned transcripts  for each cluster.
clustering was performed with pearson correlation as distance metric.
an accordingly sorted list of the clustered transcripts is provided in additional file 2: table s <dig>
time-dependent expression patterns of the prokaryotic pathogen
expression of virulence genes from the pathogen is tightly regulated, depending on the respective environment, and the ssrb response regulator of the ssra/ssrb two-component system represents one of the major virulence modulators that controls about 4% of the s. enterica genome including ssrb itself and genes in various spis  <cit> .
the first two sl <dig> pathogenicity islands  encode two distinct classical type iii secretion systems  <cit> .
these systems represent needle-like structures that inject effector proteins into the cytoplasm of host cells, which are necessary for invasion and subsequent proliferation of the pathogen.
expression of spi- <dig> genes is controlled by the ssra/ssrb along with the ompr/envz and phop/phoq two-component systems  <cit> .
the ssrb response regulator, however, exhibits the most direct effect on spi- <dig> expression, and is negatively controlled by ptsn-encoded eiiantr, a component of the nitrogen-metabolic phosphotransferase system that acts on the ssrb protein at the post-transcriptional level  <cit> .
in line with this, ptsn and ssrb display a reciprocal time-dependent expression .
the ptsn mrna is most abundant in non-interacting sl <dig> cells, absent during early interaction, and slightly upregulated in the subsequent course of infection.
conversely, the transcript encoding ssrb is most abundant during early interaction, but not detectable in non-interacting cells and mid-level interaction.
fold changes of transcripts expressed from one of the three plasmids  of non-interacting and early interacting, early and mid-level interacting, as well as mid-level and late interacting sl <dig> cells were subjected to k-means clustering  to identify transcripts with similar time-dependent expression patterns .
salmonella virulence plasmids are very stable and present in low copy numbers .
the generated clusters consequently reflect the temporal contribution of plasmid-encoded genes to the virulence program of sl <dig>  transcripts relevant for invasion of host cells and early adaption of endocytosed sl <dig>  for instance, are represented in cluster ii.
upregulation of some of these transcripts from mid-level to late interaction additionally suggests a function of the encoded proteins in later stages of infection, and amongst others this cluster includes the mrna for plasmid-encoded fimbriae d .
cluster iii, on the other hand, comprises transcripts that are especially involved in mid-level interaction.
in this cluster, the most upregulated transcript from early to mid-level interaction represents a polycistronic mrna encoding the proteins repa, tap, and repa <dig> that are required for plasmid replication.
transcripts less involved in early or mid-level interaction, but becoming successively more important with proceeding invasion and infection of host cells are represented in cluster iv.
amongst others, this cluster includes the mrna encoding plasmid sos inhibition protein a , and the most downregulated and subsequently upregulated mrna, which encodes a putative transposase .
with  <dig> out of  <dig>  cluster i comprises almost half of the identified plasmid-encoded transcripts.
however, cluster i also displays the least consistent gene expression pattern.
only four of the  <dig> transcripts are also found to be expressed in interacting cells indicating that those transcripts, which could not reliably be quantified due to insufficient sequencing depth, were grouped into this cluster.
figure 3c displays the time-dependent expression profiles of the mrnas encoding repa, tap, and repa <dig> as well as pefd in more detail.
the fimbrial chaperone protein pefd is involved in the f1–g <dig> short -assisted assembly of thick rigid mono-adhesive pili, which represent adhesive organelles necessary for bacterial attachment to target cells  <cit> .
in line with this, pefd expression is not detectable in non-interacting cells, highly induced in early interaction, and non-detectable during mid-level interaction again.
the transcript encoding repa, tap, and repa <dig>  on the other hand, is expressed in non-interacting cells, completely absent during early interaction, highly abundant in the course of intracellular replication during mid-level interaction, and finally less abundant in late interaction again.
the identified expression patterns of interacting pathogen cells reflect previous reports of infection-related gene expression and function, which corroborates the potential of dual 3’seq to reliably assess the transcriptome of prokaryotic cells during interaction, including quantification of polycistronic transcripts.
host cell responses to bacterial infection
innate immune recognition of microbial components relies on germline-encoded pattern-recognition receptors  that recognize pathogen-associated molecular patterns  of foreign cells such as bacterial lipids, lipoproteins, proteins, nucleic acids and lipopolysaccharides , the constituents of the outer membrane of gram-negative bacteria  <cit> .
toll-like receptors  represent the most important family of membrane-bound prrs, and among the ten functional tlrs in humans tlr <dig> is unique in its ability to induce two distinct signaling pathways controlled by the tirap-myd <dig> and tram-trif pairs of adaptor proteins  <cit> .
binding of bacterial lps to md <dig> and tlr <dig> induces tirap-myd88-dependent signaling at the plasma membrane, and subsequently activates the tram-trif pathway after endocytosis of tlr <dig>  in hela cells, however, the genes encoding md <dig> as well as tlr <dig> are not transcribed  <cit>  and consequently immune recognition of the pathogen must involve other prrs.
the mrnas encoding several members of the rel/nf-κb transcription factor and nf-κb inhibitor family display a marked response in gene expression upon infection .
especially, the differential expression patterns of mrnas encoding inhibitor proteins  point to an increase in iκbs marked for degradation subsequent to activation of the pathway with the onset of infection.
except for a general trend to downregulation from mid-level to late interaction, expression levels of transcripts coding for other members of the tlr <dig> signaling cascade are largely unaffected in the different interaction stages, which is in line with the fact that most of the encoded proteins are not degraded or depleted from a cytosolic pool in the context of signal transmission.
besides tlr <dig> and tlr7- <dig> that are exclusively expressed in intracellular vesicles, host cells harbor additional classes of cytosolic prrs, including rig-i-like receptors  and nod-like receptors   <cit> .
in contrast to tlrs, rlrs are present in the cytosol of all cell types, and activation by viral, but also bacterial nucleic acids results in induced expression of type i interferon and cytokines  <cit> .
the transcript encoding lgp <dig>  is strongly upregulated in the first and last infection stage, while its expression is downregulated from early to mid-level interaction .
the differential time-dependent expression of this transcript, in contrast to other rlr-encoding mrnas as mda <dig>  and rig-i , suggests a pathogen-induced activation of lgp <dig> and points to a functional role of the protein in early and late responses to infection.
the nlr family comprises more than  <dig> members for recognition of various pamps, and the identified expression profiles of corresponding mrnas differ accordingly.
the transcript encoding naip displays the strongest upregulation from non-interacting to early interacting cells, and is additionally upregulated from mid-level to late interaction.
time-dependent expression of this mrna consequently resembles the pattern of dhx <dig>  again suggesting a functional implication of the encoded protein in response to infection.figure  <dig> differential expression of membrane-bound and cytosolic prrs in response to infection.
time-dependent log <dig> fold changes of mrnas from non-interacting and interacting cells are depicted for the listed gene products, and interactions between the encoded proteins are indicated by arrows.
differentially expressed mrnas encoding proteins involved in early  and late  nf-κb signaling.
binding of bacterial lps to tlr <dig> activates two distinct signaling cascades that lead to expression of inflammatory cytokines via activation of the rel/nf-κb transcription factor family.
activation of myd <dig> results in signal transmission via interleukin- <dig> receptor-associated kinases  that in turn activate traf <dig>  this part of the signaling pathway, however, is not functional in hela cells  because of lacking md <dig> expression.
traf <dig> signaling results in activation of the kinase tak <dig> complex that phosphorylates nemo , which in turn leads to phosphorylation of the inhibitor subunits of nf-κb  via canonical ikk-β.
phosphorylated iκbs are degraded by the ubiquitin-dependent pathway, thereby releasing nf-κb, which subsequently translocates to the nucleus to function as a transcription factor .
trif-dependent signaling is also stimulated via activation of tlr <dig> by binding of double stranded rna, and enhances activation of the kinase tak <dig> complex by formation of a multimeric protein signaling complex comprising traf <dig>  tradd, pellino- <dig>  and rip <dig> .
differential expression of other cytosolic prrs that modulate the nf-κb response by cross-signaling.
nod-like receptors comprise four subfamilies, and expression profiles of the mrnas encoding naip , nod <dig>  and nlrc <dig>  as well as nlrp <dig> and nlrp <dig>  are shown along with those from the rig-i-like receptors lgp <dig> , rig-i , and mda <dig> .
distinct library preparation of the poly+ and poly− fraction by dual 3’seq allows for discriminating polyadenylated  from non-polyadenylated  microrna  precursors.
a previous study of the sl1344-induced microrna response of hela cells identified the let- <dig> family as an important regulator of major cytokines  <cit> .
in line with this, the host gene encoding mirna let-7 g  is expressed strongest in non-interacting cells and continuously repressed throughout the interaction stages .
while the level of let-7 g pre-mirna is relatively stable during interaction, the pri-mirna abundance is halved with each point of time post infection.
mir- <dig> is involved in the response of human cells to infection with leishmania major  <cit> , and the fact that the pri-mirna encoded by mir210hg is most abundant during late interaction of human host and sl <dig> cells  suggests an additional role of this mirna in the context of bacterial infection.
taken together, deepsupersage-coupled dual 3’seq provides insights into host cell responses to bacterial infection with minimal sequencing efforts not only for the protein-coding transcriptome, but also for ncrnas such as mirna precursors.
validation of time-dependent expression patterns and assessment of the quantification accuracy of dual 3’seq by quantitative real-time pcr
the identified time-dependent expression patterns from pathogen and host cells were confirmed by miqe-conform qrt-pcr  <cit>  with three independent biological replicates.
gapd, b2m and rpl13a served as reference genes for normalization of the transcript abundances in host cells.
out of seven previously determined candidate genes for normalization of prokaryotic transcript abundances nagd, ndh, rpod, and trma exhibited the highest expression stabilities between pools of reverse-transcribed cdna from the three biological replicates of each interaction stage as well as non-interacting cells and between the three distinct biological replicates of non-interacting cells .table  <dig>
relative expression stabilities of seven candidate reference genes from prokaryotic cells according to bestkeeper, normfinder, and genorm
evaluation method
candidate reference gene stability*
trma
,
nagd
,
ndh
,
rpod
, ygfe, rpso, rpst
rpod
,
ndh
,
nagd
,
trma
, ygfe, rpso, rpst
trma
/
nagd
,
ndh
,
rpod
, ygfe, rpso, rpst
*gene stability decreasing from left to right; the employed reference genes for quantification of the time-dependent expression from prokaryotic transcripts are indicated by bold letters.
please consult the text for more information.
the biological variance of  <dig> target mrnas involved in host cell responses to bacterial infection is relatively small compared to the variance of  <dig> target transcripts from the pathogen .
an exact qrt-pcr-based quantification of prokaryotic transcripts is apparently impeded by the high ct values that are associated with the relatively low ratio of pathogen-derived transcripts in all interaction stages.
clustering of the time-dependent differential expression from the  <dig> human mrnas either determined by dual 3’seq coupled to deepsupersage or by qrt-pcr, on the other hand, reveals a good correlation especially for the mrnas encoded by nfkb <dig>  nfkb <dig>  nfkbie and relb .
according to this, the initial immune response of host cells involves activation of p <dig>  and p <dig>  that act as transcriptional activators or repressors depending on the respective dimerization partner, along with i-rel , which exhibits the most prominent differential expression of all  <dig> mrnas.
since a consistent differential expression of these mrnas across all three biological replicates is only present from early to mid-level and from mid-level to late interaction, the initial host cell response to sl <dig> infection seems to occur in between half an hour and four hours post infection.
principal component analysis of the  <dig> human mrnas demonstrates that the time-dependent differential expression measured by qrt-pcr of the biological replicates is similar to the respective expression patterns determined by dual 3’seq , and allows for discriminating the different interaction stages from each other .figure  <dig> confirmation of time-dependent expression patterns from pathogen as well as host cells and evaluation of dual 3’seq quantification accuracy by quantitative real-time pcr.
biological variance of the targeted transcripts across three independent biological replicates.
means with sd of the ∆ct values from non-interacting cells  and cells from early , mid-level  and late interaction  are shown for  <dig> mrnas involved in host cell responsaes to bacterial infection  as well as  <dig> prokaryotic transcripts .
additional coding sequences of polycistronic transcripts from the prokaryote are indicated by brackets.
∆ct values are plotted on an inverted ordinate, where a low figure reflects an increased abundance of the respective transcript and vice versa.
for each graph, the lowest ct value was arbitrarily shifted to zero followed by an according normalization of the other values.
heat map of time-dependent log <dig> fold changes and corresponding ∆∆ct values from early and mid-level interacting as well as mid-level and late interacting host cells.
the individual expression ratios determined by qrt-pcr from three biological replicates  together with the corresponding ratios identified by deepsupersage  are depicted after hierarchical clustering of all samples based on uncentered pearson correlation.
upregulated transcripts are represented in red, while downregulated transcripts are shown in green.
principal component analysis of the data included in the heat map along with the corresponding log <dig> fold changes from non-interacting and early interacting cells .
scatterplot depicting the relative quantification accuracy of dual 3’seq in comparison to probe-based qrt-pcr.
the expression ratios of  <dig> polyadenylated mrnas and their non-polyadenylated degradation intermediates in early interacting cells determined either by deepsupersage  or mace  is plotted against the respective qrt-pcr-based ratios .
an optimal correlation is indicated by the diagonal line.
quantification accuracy of dual 3’seq coupled either to deepsupersage or alternatively to mace was assessed by probe-based qrt-pcr of  <dig> mrnas with dnase-treated and rrna-depleted total rna from library preparation of the early interaction stage.
the ratios between functional mrnas and their non-polyadenylated degradation intermediates were determined in one-step qrt-pcr reactions and plotted against the identified ratios according to deepsupersage or mace, respectively .
even though two of the targets  are not present in the sequencing data of deepsupersage-coupled library preparation, the pearson correlation is better than the correlation of mace-based dual 3’seq .
four out of the  <dig> targets captured by deepsupersage exhibit an almost optimal correlation with qrt-pcr, and especially the transcript encoded by slc7a <dig> is found to be more abundant in the non-polyadenylated transcriptome by both methods.
CONCLUSIONS
the identified ratios between human transcripts and transcripts originating from sl <dig> in interacting cells demand extensive sequencing to sufficiently cover the transcriptome of interacting pathogen cells without their prior enrichment.
compared to the eukaryote, reads from sl <dig> exhibit a maximal ratio of  <dig> to  <dig>  additionally, sl1344-derived reads in the poly+ fraction are notably less abundant  compared to the corresponding poly− fraction .
the latter fraction is crucial for functional analysis of differentially expressed transcripts that encode proteins involved in infection, since polyadenylated transcripts from the prokaryote merely represent degradation intermediates.
distinct library preparation with the poly+ and poly− fraction of interacting cells via dual 3’seq preserves the relatively high percentages of sl1344-derived transcripts in the non-polyadenylated fraction of host cells, and consequently permits a more detailed analysis of the functional, non-polyadenylated prokaryotic transcriptome during interaction.
in contrast, combined library preparation of both fractions inevitably leads to a reduced ratio of pathogen transcripts due to the immense diversity of the polyadenylated host cell transcriptome on the one, and the extremely small number of polyadenylated intermediates from the pathogen, on the other hand.
reciprocal expression of ptsn and ssrb proves that polycistronic transcripts from the prokaryote were accurately quantified, since ptsn is transcribed as polycistronic mrna that encodes ten other gene products, while ssrb is expressed as a monocistronic transcript.
the efficiency of in-silico separation of the two interacting organisms can be assessed through the number of false positive annotated reads in the libraries of non-interacting cells .
although, the poly+ fraction of non-interacting sl <dig> cells contains a relatively high percentage of incorrectly aligned reads , the corresponding poly− fraction only comprises  <dig> % of false positive reads.
with  <dig> % in the poly+ and  <dig> % in the poly− fraction, libraries from non-interacting host cells comprise even less incorrectly aligned reads.
thus, the high ratio of false positive reads in the poly+ fraction of non-interacting sl <dig> cells is not representative for the general efficiency of the employed annotation procedure, but rather caused by the extremely low coverage of this library and sequencing of degraded intermediates.
mapping of 3’ reads with a length of more than  <dig> nucleotides did not substantially improve the ratio of transcripts that aligned to both the human and the sl <dig> transcriptome as can be inferred from the numbers of reads excluded during multi-step annotation .
the corresponding percentages of excluded reads from early interaction libraries prepared either with deepsupersage or mace are relatively similar + and  <dig> % vs.  <dig> % in the poly− fraction, respectively).
however, generation of tags via deepsupersage depends on the presence of a recognition site for the anchoring enzyme, while mace generates a tag out of every transcript.
the comparison of both methods with qrt-pcr confirms that not all expressed transcripts are represented in the deepsupersage dataset, but also that the captured transcripts are quantified accurately by deepsupersage-coupled dual 3’seq.
transcriptome profiling by deepsupersage consequently provides a reliable representation of the differential expression from interacting pro- and eukaryotic cells, but it also involves a much more elaborate library preparation procedure compared to mace.
another advantage of the latter technique is its potential to capture the expression of small ncrnas, since it does not necessarily involve size-selection .
additionally, mace can be adjusted to any desired sequencing length, which allows for identification of alternative polyadenylation sites in eukaryotes  <cit> .
for functional analysis, the trade-off between sequencing depth and the number of biological replicates for each condition must be carefully balanced  <cit> .
consequently, every decrease in the necessary amount of sequencing per sample allows for a more detailed time-dependent analysis and more biological replicates per condition.
dual 3’seq employing deepsupersage significantly reduces the complexity of the involved transcriptomes.
even though the prokaryotic transcriptome of interacting cells was not completely covered, our approach provides first insights into the pathogenicity-related gene expression of sl <dig> and corresponding host cell responses.
the time-dependent expression patterns of spi and plasmid-encoded transcripts in interacting sl <dig> cells reflect a coordinated activation of virulence genes in the course of infection, which in turn elicits corresponding host responses in transcription as exemplarily shown for the signaling cascades leading to expression of inflammatory cytokines.

 interactions between eu- and prokaryotic cells are frequent, multifaceted events ranging from symbiotic synergy such as symbiotic nitrogen fixation in legumes or fermentation by gastrointestinal bacteria to pathogenic interference, for instance, in the course of salmonellosis. this interplay of organisms requires mutual signaling mechanisms and a continuous adaptation of the metabolism of the involved cells to varying environmental conditions. consequently, programmed expression patterns have to be induced to continuously readjust the proteome and metabolome of both cell types. the characterization of corresponding time-dependent expression patterns allows for a deeper understanding of the underlying molecular processes, and was the focus of numerous studies, but until recently gene expression profiling emphasized either the host cell or the prokaryotic transcriptome. salmonella represents a genus of gram-negative and facultative anaerobic enterobacteria and is closely related to the genus escherichia. the human-virulent pathogen is a model organism for characterization of host-pathogen interactions, and is associated with a variety of diseases including gastroenteritis and enteric fever. one of the first attempts to profile genome-wide expression changes in host cells during interaction of pro- and eukaryotes was carried out with s. enterica and human intestinal epithelial cells. conversely, the transcriptome of s. enterica became subject to several studies of host-pathogen interactions after completion of the genome sequences of s. enterica serotype typhi ct and serotype typhimurium lt, which complemented some of the previously determined host responses. in the meantime, next-generation sequencing -coupled transcription profiling techniques emerged as the principal tools to interrogate gene expression, and especially whole transcriptome shotgun sequencing has considerably contributed to our understanding of prokaryotic transcriptomes. nonetheless, simultaneous transcription profiling without prior disruption of the interaction remains technically challenging, and thus characterization of disease-related expression patterns in interacting eu- and prokaryotic cells is inevitably linked to comprehensive sequencing efforts. here we present dual 3’seq, a tag-based, ngs-coupled method that allows for simultaneous transcription profiling of interacting pro- and eukaryotes without physical separation of the interacting cells. compared to rna-seq, the reduction in complexity of tag-based approaches significantly decreases the required sequencing depth for a good coverage of both the pro- and eukaryotic transcriptomes, which is a prerequisite for profiling of low abundant pathogen-derived transcripts. additionally, only a single tag is generated out of each transcript, which facilitates unequivocal quantification of reads from a specific rna without sacrificing qualitative information of pathogen-derived transcripts, since prokaryotes lack alternative splicing events. deepsupersage and mace represent two established ngs-coupled transcriptome profiling techniques that generate exactly one tag out of the 3′ end of every transcript. while deepsupersage yields a nucleotide tag that is specifically located within the 3′ end depending on the presence of an according restriction site for the employed anchoring enzyme, mace generates randomly distributed tags out of the last hundreds of bases. mace consequently allows for preparation of libraries with varying read lengths and provides additional information regarding transcripts that do not possess an according restriction site for anchoring. in order to assess the respective efficiencies in transcriptome profiling of cultivated and interacting s. enterica typhimurium and human host cells, we combined the dual 3’seq approach with both protocols. human epithelial cells were infected with the invasive pathogen s. enterica typhimurium sl, and interacting cells were screened for differentially expressed transcripts at several points of time post infection to provide an overview of the transcriptional processes during invasion of the intestinal epithelium as one of the first steps in emerging salmonellosis. the combination of the published sl transcriptome with the operon structure identified by differential rna-seq allowed for accurate quantification of polycistronically transcribed genes from the prokaryote, and the corresponding expression profiles provide a basis for time-dependent analysis of disease-related transcripts, despite the fact that transcripts from the pathogen are of extremely low abundance during interaction.  dual 3’seq of interacting pro- and eukaryotes in order to reduce the sequencing depth for simultaneous characterization of the transcriptomes of interacting pro- and eukaryotic cells without prior separation, we established the new transcriptome profiling technique dual 3’seq. we employed this technique to identify time-dependent transcription patterns of interacting hela-s and sl cells subsequent to infection of the host cells. total rna from cultured and interacting cells was isolated after different points of time to interrogate the molecular events during early, mid-level and late interaction stages of infection. the first interaction point was chosen as early as technically possible to capture the initial pathogen responses after invasion of the host cells, while mid-level and late interaction allow for investigation of the early and late host cell response as well as intracellular pathogen replication.figure scheme of dual 3’seq library preparation and bioinformatic processing of the generated sequencing data. total rna was size-selected subsequent to dnase i digestion of remaining dna in the isolate. following rrna depletion, the rna was split into the poly+ and poly− fraction by oligo capture to separate the polyadenylated and functional mrnas of eukaryotic cells from the non-polyadenylated transcripts that represent the functional transcriptome of prokaryotes. ensuing in-vitro polyadenylation of the poly− fraction, both fractions were subjected to oligo-based reverse transcription. the generated cdna was fragmented according to two established 3′ transcriptome profiling techniques. deepsupersage tags were generated via cleavage of rnas by the anchoring enzyme nlaiii and subsequent digestion using ecop15i, while mace involved random fragmentation for generation of tags. 3′ fragments were enriched by binding to a streptavidin matrix and ligated to a sequencing adaptor. adaptor-ligated fragments were pcr-amplified using genxpro’s truequant technology for pcr-bias free amplification, page-purified, and sequenced on the illumina hiseq platform. barcoded reads were allocated to their respective library, filtered for pcr-derived reads, and trimmed for high-quality sequences. afterwards, reads were annotated to a combined reference comprising the transcriptome and genome sequences of sl and human host cells in a multi-step procedure. reads uniquely mapped to one of both organisms were combined to three distinct expression matrices for functional analysis of the poly− transcriptome from pathogen and host cells as well as the poly+ fraction of the host cells. for each expression matrix, annotated reads were quantified and median-normalized using deseq, followed by pair-wise, time-dependent comparison of the different interaction stages. statistical significance was subsequently corrected for multiple testing according to benjamini and hochberg. the generated sequencing data from the poly− fractions comprise between and million high-quality and pcr duplicate-free reads, while the number of reads in the corresponding poly+ libraries ranges from almost million down to ~ reads in the sl library. the relatively low number of filtered poly+ reads from cultivated and especially interacting sl cells reflects post-transcriptional regulation in prokaryotes, where polyadenylation facilitates degradation of rna. additionally, the low abundance of polyadenylated intermediates from the prokaryote is corroborated by the low ratios of sl1344-derived transcripts in all interaction stages, which impedes a functional analysis of this fraction. prokaryotic reads in libraries of the poly− fraction are about ten times more abundant. compared to the ratios of the corresponding poly+ libraries, non-polyadenylated sl reads are about 160-fold, 50-fold and 10-fold higher from early to late interaction, respectively. this increase in polyadenylated reads during interaction suggests an increased turnover of transcripts probably caused by the growing number of pathogen cells in the confined environment of hela cell culture. discarding ambiguously annotated poly− reads that aligned to more than one gene, or 89% of the known loci in sl are transcribed in cultivated or interacting cells, while or 75% of protein-coding loci are represented in the polyadenylated transcriptomes of human host cells. over 95% of the transcripts in human poly+ libraries encode proteins, and mrnas encoding ribosomal proteins account for a large proportion of these. with an mrna content ranging from 13% to 20% protein-coding transcripts are much less abundant in the corresponding poly− libraries, which mostly comprise small nucleolar rnas, small nuclear rnas and other non-coding rnas.figure annotation statistics for the generated expression matrices and differentially expressed genes of interacting pathogen and host cells. read numbers annotated to human and pathogen cells or both for all sequenced libraries, respectively. the bar graph shows log10-transformed numbers, while the original number of reads is depicted in a circle graph. classification of unambiguously annotated transcripts from host cells into mrnas, non-coding rnas, bivalent transcripts encoding mrnas and non-coding rnas, long non-coding rnas, histone-encoding mrnas, and other transcript classes for the poly+ and poly− transcriptome, respectively. mrnas encoding cytoplasmic or mitochondrial ribosomal proteins are depicted individually for the poly+ and combined for the poly− transcriptome that additionally comprises small nucleolar rnas and small nuclear rnas. respective abundances are shown for non-interacting as well as interacting cells. the different interaction stages prepared with 3’seq coupled to deepsupersage are depicted from inside outwards along with the early interaction library that was prepared using mace. venn diagram showing the number of commonly expressed genes for all combinations of libraries corresponding to one of the expression matrices. ambiguously annotated reads are not included. pair-wise comparison of differentially expressed transcripts within the expression matrices along with the number of exclusively expressed transcripts for each library. listed are unambiguously annotated reads with a log fold change stronger than | | and an fdr-corrected p-value below. up- and downregulated transcripts are indicated by corresponding arrows for the library on the ordinate, respectively. differential gene expression during interaction the intersections of unambiguously annotated rnas from the poly− transcriptome of sl cells as well as both transcriptome fractions of human host cells reveal distinct overlaps of the expressed genes. the four sl libraries comprise ~ commonly expressed genes, while non-interacting and late interacting cells share the highest number with more than. unsurprisingly, cells from early and mid-level interaction express fewer genes in common with non-interacting cells than those from the late interaction stage when infection is well established. conversely, the number of commonly expressed genes in the poly+ transcriptome of host cells is highest between early and mid-level interaction. with  ±  genes, the numbers of overlapping genes in the poly− fraction of host cells spread notably less than in the poly+ fraction, which suggests a more diverse effect of infection on the polyadenylated transcriptome of host cells. although the poly− fraction of host cells generally comprises more significantly up and downregulated transcripts than the poly+ fraction, expression of polyadenylated transcripts during infection is more affected, since mrna degradation naturally leads to corresponding changes in the non-polyadenylated transcriptome. in fact, more than 80% of the non-polyadenylated, significantly differentially expressed transcripts in all the pairwise comparisons represent protein-coding mrnas, even though these rnas only constitute up to a fifth of their respective poly− fractions in total. expression of stm and stm is non-detectable except for the transcriptome of early interacting sl cells, while the transcripts encoded by gst, ydin and stm are only present in mid-level interaction. stm is part of salmonella pathogenicity island, and substantially contributes to bacterial fitness by promoting transcription of spi- together with ssrb. stm and stm are both encoded within gifsy- a lambdoid prophage that promotes virulence in addition to the remnant phage encoded by spi-. the time-dependent differential expression of all sequenced genes from the prokaryote reveals characteristic patterns for each interaction stage. the expression of many spi-encoded and phage-derived transcripts from the chromosome as well as plasmid-derived transcripts is induced or repressed in the course of infection, and analysis of these disease-related expression patterns allows for dissolving the time-dependent components of sl gene expression during different stages of infection.figure differential expression of interacting pathogen cells. circos plot of sl gene expression for chromosomally-encoded genes as well as plasmid-derived transcripts. the respective coordinates are given in megabases along with the corresponding log fold changes of non-interacting and early interacting, early and mid-level interacting, as well as mid-level and late interacting cells. salmonella pathogenicity islands are highlighted in orange, while other important loci are shown in green. detailed time-dependent expression profiles of selected sl genes. the median-normalized transcript abundance is shown for non-interacting cells along with those from the different interaction stages. ptsn is transcribed as polycistronic mrna from the chromosome that comprises ten other genes. the mrna encoding pefd and the polycistronic mrna coding for repa, tap, and repa are both transcribed from intracellular plasmid psltsl k-means clustering of the time-dependent fold changes from plasmid-encoded transcripts. the corresponding log fold changes were assigned to one of four clusters based on their time-dependent expression patterns. the closest centroid is shown along with the assigned transcripts for each cluster. clustering was performed with pearson correlation as distance metric. an accordingly sorted list of the clustered transcripts is provided in additional file 2: table s time-dependent expression patterns of the prokaryotic pathogen expression of virulence genes from the pathogen is tightly regulated, depending on the respective environment, and the ssrb response regulator of the ssra/ssrb two-component system represents one of the major virulence modulators that controls about 4% of the s. enterica genome including ssrb itself and genes in various spis. the first two sl pathogenicity islands encode two distinct classical type iii secretion systems. these systems represent needle-like structures that inject effector proteins into the cytoplasm of host cells, which are necessary for invasion and subsequent proliferation of the pathogen. expression of spi- genes is controlled by the ssra/ssrb along with the ompr/envz and phop/phoq two-component systems. the ssrb response regulator, however, exhibits the most direct effect on spi- expression, and is negatively controlled by ptsn-encoded eiiantr, a component of the nitrogen-metabolic phosphotransferase system that acts on the ssrb protein at the post-transcriptional level. in line with this, ptsn and ssrb display a reciprocal time-dependent expression. the ptsn mrna is most abundant in non-interacting sl cells, absent during early interaction, and slightly upregulated in the subsequent course of infection. conversely, the transcript encoding ssrb is most abundant during early interaction, but not detectable in non-interacting cells and mid-level interaction. fold changes of transcripts expressed from one of the three plasmids of non-interacting and early interacting, early and mid-level interacting, as well as mid-level and late interacting sl cells were subjected to k-means clustering to identify transcripts with similar time-dependent expression patterns. salmonella virulence plasmids are very stable and present in low copy numbers. the generated clusters consequently reflect the temporal contribution of plasmid-encoded genes to the virulence program of sl transcripts relevant for invasion of host cells and early adaption of endocytosed sl for instance, are represented in cluster ii. upregulation of some of these transcripts from mid-level to late interaction additionally suggests a function of the encoded proteins in later stages of infection, and amongst others this cluster includes the mrna for plasmid-encoded fimbriae d. cluster iii, on the other hand, comprises transcripts that are especially involved in mid-level interaction. in this cluster, the most upregulated transcript from early to mid-level interaction represents a polycistronic mrna encoding the proteins repa, tap, and repa that are required for plasmid replication. transcripts less involved in early or mid-level interaction, but becoming successively more important with proceeding invasion and infection of host cells are represented in cluster iv. amongst others, this cluster includes the mrna encoding plasmid sos inhibition protein a, and the most downregulated and subsequently upregulated mrna, which encodes a putative transposase. with out of cluster i comprises almost half of the identified plasmid-encoded transcripts. however, cluster i also displays the least consistent gene expression pattern. only four of the transcripts are also found to be expressed in interacting cells indicating that those transcripts, which could not reliably be quantified due to insufficient sequencing depth, were grouped into this cluster. figure 3c displays the time-dependent expression profiles of the mrnas encoding repa, tap, and repa as well as pefd in more detail. the fimbrial chaperone protein pefd is involved in the f1–g short -assisted assembly of thick rigid mono-adhesive pili, which represent adhesive organelles necessary for bacterial attachment to target cells. in line with this, pefd expression is not detectable in non-interacting cells, highly induced in early interaction, and non-detectable during mid-level interaction again. the transcript encoding repa, tap, and repa on the other hand, is expressed in non-interacting cells, completely absent during early interaction, highly abundant in the course of intracellular replication during mid-level interaction, and finally less abundant in late interaction again. the identified expression patterns of interacting pathogen cells reflect previous reports of infection-related gene expression and function, which corroborates the potential of dual 3’seq to reliably assess the transcriptome of prokaryotic cells during interaction, including quantification of polycistronic transcripts. host cell responses to bacterial infection innate immune recognition of microbial components relies on germline-encoded pattern-recognition receptors that recognize pathogen-associated molecular patterns of foreign cells such as bacterial lipids, lipoproteins, proteins, nucleic acids and lipopolysaccharides, the constituents of the outer membrane of gram-negative bacteria. toll-like receptors represent the most important family of membrane-bound prrs, and among the ten functional tlrs in humans tlr is unique in its ability to induce two distinct signaling pathways controlled by the tirap-myd and tram-trif pairs of adaptor proteins. binding of bacterial lps to md and tlr induces tirap-myd88-dependent signaling at the plasma membrane, and subsequently activates the tram-trif pathway after endocytosis of tlr in hela cells, however, the genes encoding md as well as tlr are not transcribed and consequently immune recognition of the pathogen must involve other prrs. the mrnas encoding several members of the rel/nf-κb transcription factor and nf-κb inhibitor family display a marked response in gene expression upon infection. especially, the differential expression patterns of mrnas encoding inhibitor proteins point to an increase in iκbs marked for degradation subsequent to activation of the pathway with the onset of infection. except for a general trend to downregulation from mid-level to late interaction, expression levels of transcripts coding for other members of the tlr signaling cascade are largely unaffected in the different interaction stages, which is in line with the fact that most of the encoded proteins are not degraded or depleted from a cytosolic pool in the context of signal transmission. besides tlr and tlr7- that are exclusively expressed in intracellular vesicles, host cells harbor additional classes of cytosolic prrs, including rig-i-like receptors and nod-like receptors. in contrast to tlrs, rlrs are present in the cytosol of all cell types, and activation by viral, but also bacterial nucleic acids results in induced expression of type i interferon and cytokines. the transcript encoding lgp is strongly upregulated in the first and last infection stage, while its expression is downregulated from early to mid-level interaction. the differential time-dependent expression of this transcript, in contrast to other rlr-encoding mrnas as mda and rig-i, suggests a pathogen-induced activation of lgp and points to a functional role of the protein in early and late responses to infection. the nlr family comprises more than members for recognition of various pamps, and the identified expression profiles of corresponding mrnas differ accordingly. the transcript encoding naip displays the strongest upregulation from non-interacting to early interacting cells, and is additionally upregulated from mid-level to late interaction. time-dependent expression of this mrna consequently resembles the pattern of dhx again suggesting a functional implication of the encoded protein in response to infection.figure differential expression of membrane-bound and cytosolic prrs in response to infection. time-dependent log fold changes of mrnas from non-interacting and interacting cells are depicted for the listed gene products, and interactions between the encoded proteins are indicated by arrows. differentially expressed mrnas encoding proteins involved in early and late nf-κb signaling. binding of bacterial lps to tlr activates two distinct signaling cascades that lead to expression of inflammatory cytokines via activation of the rel/nf-κb transcription factor family. activation of myd results in signal transmission via interleukin- receptor-associated kinases that in turn activate traf this part of the signaling pathway, however, is not functional in hela cells because of lacking md expression. traf signaling results in activation of the kinase tak complex that phosphorylates nemo, which in turn leads to phosphorylation of the inhibitor subunits of nf-κb via canonical ikk-β. phosphorylated iκbs are degraded by the ubiquitin-dependent pathway, thereby releasing nf-κb, which subsequently translocates to the nucleus to function as a transcription factor. trif-dependent signaling is also stimulated via activation of tlr by binding of double stranded rna, and enhances activation of the kinase tak complex by formation of a multimeric protein signaling complex comprising traf tradd, pellino- and rip. differential expression of other cytosolic prrs that modulate the nf-κb response by cross-signaling. nod-like receptors comprise four subfamilies, and expression profiles of the mrnas encoding naip, nod and nlrc as well as nlrp and nlrp are shown along with those from the rig-i-like receptors lgp, rig-i, and mda. distinct library preparation of the poly+ and poly− fraction by dual 3’seq allows for discriminating polyadenylated from non-polyadenylated microrna precursors. a previous study of the sl1344-induced microrna response of hela cells identified the let- family as an important regulator of major cytokines. in line with this, the host gene encoding mirna let-7 g is expressed strongest in non-interacting cells and continuously repressed throughout the interaction stages. while the level of let-7 g pre-mirna is relatively stable during interaction, the pri-mirna abundance is halved with each point of time post infection. mir- is involved in the response of human cells to infection with leishmania major, and the fact that the pri-mirna encoded by mir210hg is most abundant during late interaction of human host and sl cells suggests an additional role of this mirna in the context of bacterial infection. taken together, deepsupersage-coupled dual 3’seq provides insights into host cell responses to bacterial infection with minimal sequencing efforts not only for the protein-coding transcriptome, but also for ncrnas such as mirna precursors. validation of time-dependent expression patterns and assessment of the quantification accuracy of dual 3’seq by quantitative real-time pcr the identified time-dependent expression patterns from pathogen and host cells were confirmed by miqe-conform qrt-pcr with three independent biological replicates. gapd, b2m and rpl13a served as reference genes for normalization of the transcript abundances in host cells. out of seven previously determined candidate genes for normalization of prokaryotic transcript abundances nagd, ndh, rpod, and trma exhibited the highest expression stabilities between pools of reverse-transcribed cdna from the three biological replicates of each interaction stage as well as non-interacting cells and between the three distinct biological replicates of non-interacting cells.table relative expression stabilities of seven candidate reference genes from prokaryotic cells according to bestkeeper, normfinder, and genorm evaluation method candidate reference gene stability* trma, nagd, ndh, rpod, ygfe, rpso, rpst rpod, ndh, nagd, trma, ygfe, rpso, rpst trma / nagd, ndh, rpod, ygfe, rpso, rpst *gene stability decreasing from left to right; the employed reference genes for quantification of the time-dependent expression from prokaryotic transcripts are indicated by bold letters. please consult the text for more information. the biological variance of target mrnas involved in host cell responses to bacterial infection is relatively small compared to the variance of target transcripts from the pathogen. an exact qrt-pcr-based quantification of prokaryotic transcripts is apparently impeded by the high ct values that are associated with the relatively low ratio of pathogen-derived transcripts in all interaction stages. clustering of the time-dependent differential expression from the human mrnas either determined by dual 3’seq coupled to deepsupersage or by qrt-pcr, on the other hand, reveals a good correlation especially for the mrnas encoded by nfkb nfkb nfkbie and relb. according to this, the initial immune response of host cells involves activation of p and p that act as transcriptional activators or repressors depending on the respective dimerization partner, along with i-rel, which exhibits the most prominent differential expression of all mrnas. since a consistent differential expression of these mrnas across all three biological replicates is only present from early to mid-level and from mid-level to late interaction, the initial host cell response to sl infection seems to occur in between half an hour and four hours post infection. principal component analysis of the human mrnas demonstrates that the time-dependent differential expression measured by qrt-pcr of the biological replicates is similar to the respective expression patterns determined by dual 3’seq, and allows for discriminating the different interaction stages from each other.figure confirmation of time-dependent expression patterns from pathogen as well as host cells and evaluation of dual 3’seq quantification accuracy by quantitative real-time pcr. biological variance of the targeted transcripts across three independent biological replicates. means with sd of the ∆ct values from non-interacting cells and cells from early, mid-level and late interaction are shown for mrnas involved in host cell responsaes to bacterial infection as well as prokaryotic transcripts. additional coding sequences of polycistronic transcripts from the prokaryote are indicated by brackets. ∆ct values are plotted on an inverted ordinate, where a low figure reflects an increased abundance of the respective transcript and vice versa. for each graph, the lowest ct value was arbitrarily shifted to zero followed by an according normalization of the other values. heat map of time-dependent log fold changes and corresponding ∆∆ct values from early and mid-level interacting as well as mid-level and late interacting host cells. the individual expression ratios determined by qrt-pcr from three biological replicates together with the corresponding ratios identified by deepsupersage are depicted after hierarchical clustering of all samples based on uncentered pearson correlation. upregulated transcripts are represented in red, while downregulated transcripts are shown in green. principal component analysis of the data included in the heat map along with the corresponding log fold changes from non-interacting and early interacting cells. scatterplot depicting the relative quantification accuracy of dual 3’seq in comparison to probe-based qrt-pcr. the expression ratios of polyadenylated mrnas and their non-polyadenylated degradation intermediates in early interacting cells determined either by deepsupersage or mace is plotted against the respective qrt-pcr-based ratios. an optimal correlation is indicated by the diagonal line. quantification accuracy of dual 3’seq coupled either to deepsupersage or alternatively to mace was assessed by probe-based qrt-pcr of mrnas with dnase-treated and rrna-depleted total rna from library preparation of the early interaction stage. the ratios between functional mrnas and their non-polyadenylated degradation intermediates were determined in one-step qrt-pcr reactions and plotted against the identified ratios according to deepsupersage or mace, respectively. even though two of the targets are not present in the sequencing data of deepsupersage-coupled library preparation, the pearson correlation is better than the correlation of mace-based dual 3’seq. four out of the targets captured by deepsupersage exhibit an almost optimal correlation with qrt-pcr, and especially the transcript encoded by slc7a is found to be more abundant in the non-polyadenylated transcriptome by both methods.  the identified ratios between human transcripts and transcripts originating from sl in interacting cells demand extensive sequencing to sufficiently cover the transcriptome of interacting pathogen cells without their prior enrichment. compared to the eukaryote, reads from sl exhibit a maximal ratio of to additionally, sl1344-derived reads in the poly+ fraction are notably less abundant compared to the corresponding poly− fraction. the latter fraction is crucial for functional analysis of differentially expressed transcripts that encode proteins involved in infection, since polyadenylated transcripts from the prokaryote merely represent degradation intermediates. distinct library preparation with the poly+ and poly− fraction of interacting cells via dual 3’seq preserves the relatively high percentages of sl1344-derived transcripts in the non-polyadenylated fraction of host cells, and consequently permits a more detailed analysis of the functional, non-polyadenylated prokaryotic transcriptome during interaction. in contrast, combined library preparation of both fractions inevitably leads to a reduced ratio of pathogen transcripts due to the immense diversity of the polyadenylated host cell transcriptome on the one, and the extremely small number of polyadenylated intermediates from the pathogen, on the other hand. reciprocal expression of ptsn and ssrb proves that polycistronic transcripts from the prokaryote were accurately quantified, since ptsn is transcribed as polycistronic mrna that encodes ten other gene products, while ssrb is expressed as a monocistronic transcript. the efficiency of in-silico separation of the two interacting organisms can be assessed through the number of false positive annotated reads in the libraries of non-interacting cells. although, the poly+ fraction of non-interacting sl cells contains a relatively high percentage of incorrectly aligned reads, the corresponding poly− fraction only comprises % of false positive reads. with % in the poly+ and % in the poly− fraction, libraries from non-interacting host cells comprise even less incorrectly aligned reads. thus, the high ratio of false positive reads in the poly+ fraction of non-interacting sl cells is not representative for the general efficiency of the employed annotation procedure, but rather caused by the extremely low coverage of this library and sequencing of degraded intermediates. mapping of 3’ reads with a length of more than nucleotides did not substantially improve the ratio of transcripts that aligned to both the human and the sl transcriptome as can be inferred from the numbers of reads excluded during multi-step annotation. the corresponding percentages of excluded reads from early interaction libraries prepared either with deepsupersage or mace are relatively similar + and % vs. % in the poly− fraction, respectively). however, generation of tags via deepsupersage depends on the presence of a recognition site for the anchoring enzyme, while mace generates a tag out of every transcript. the comparison of both methods with qrt-pcr confirms that not all expressed transcripts are represented in the deepsupersage dataset, but also that the captured transcripts are quantified accurately by deepsupersage-coupled dual 3’seq. transcriptome profiling by deepsupersage consequently provides a reliable representation of the differential expression from interacting pro- and eukaryotic cells, but it also involves a much more elaborate library preparation procedure compared to mace. another advantage of the latter technique is its potential to capture the expression of small ncrnas, since it does not necessarily involve size-selection. additionally, mace can be adjusted to any desired sequencing length, which allows for identification of alternative polyadenylation sites in eukaryotes. for functional analysis, the trade-off between sequencing depth and the number of biological replicates for each condition must be carefully balanced. consequently, every decrease in the necessary amount of sequencing per sample allows for a more detailed time-dependent analysis and more biological replicates per condition. dual 3’seq employing deepsupersage significantly reduces the complexity of the involved transcriptomes. even though the prokaryotic transcriptome of interacting cells was not completely covered, our approach provides first insights into the pathogenicity-related gene expression of sl and corresponding host cell responses. the time-dependent expression patterns of spi and plasmid-encoded transcripts in interacting sl cells reflect a coordinated activation of virulence genes in the course of infection, which in turn elicits corresponding host responses in transcription as exemplarily shown for the signaling cascades leading to expression of inflammatory cytokines. 
eighty-nine percent of the known loci are found to be transcribed in prokaryotic cells prior or subsequent to infection of the host, while 75% of all protein-coding loci are represented in the polyadenylated transcriptomes of human host cells. the interaction of eukaryotic host and prokaryotic pathogen cells is linked to specific changes in the cellular proteome, and consequently to infection-related gene expression patterns of the involved cells. additionally, the identified expression patterns of both organisms were validated by qrt-pcr using three independent biological replicates, which confirmed that relb along with nfkb and nfkb are involved in the initial immune response of epithelial cells after infection with s. enterica typhimurium. human epithelial cells were infected with salmonella enterica typhimurium as a model system for invasion of the intestinal epithelium, and the transcriptional response of the infected host cells together with the differential expression of invading and intracellular pathogen cells was determined by dual 3’seq coupled with the next-generation sequencing-based transcriptome profiling technique deepsupersage. to simultaneously assess the transcriptomes of both organisms during their interaction we developed dual 3’seq, a tag-based sequencing protocol that allows for exact quantification of differentially expressed transcripts in interacting pro- and eukaryotic cells without prior fixation or physical disruption of the interaction. dual 3’seq was alternatively coupled to mace to assess the advantages and drawbacks of a library preparation procedure that allows for sequencing of longer fragments. annotation to reference transcriptomes comprising the operon structure of the employed s. enterica typhimurium strain allowed for in silico separation of the interacting cells including quantification of polycistronic rnas.
 the interaction of eukaryotic host and prokaryotic pathogen cells is linked to specific changes in the cellular proteome, and consequently to infection-related gene expression patterns of the involved cells. to simultaneously assess the transcriptomes of both organisms during their interaction we developed dual 3’seq, a tag-based sequencing protocol that allows for exact quantification of differentially expressed transcripts in interacting pro- and eukaryotic cells without prior fixation or physical disruption of the interaction.  human epithelial cells were infected with salmonella enterica typhimurium as a model system for invasion of the intestinal epithelium, and the transcriptional response of the infected host cells together with the differential expression of invading and intracellular pathogen cells was determined by dual 3’seq coupled with the next-generation sequencing-based transcriptome profiling technique deepsupersage. annotation to reference transcriptomes comprising the operon structure of the employed s. enterica typhimurium strain allowed for in silico separation of the interacting cells including quantification of polycistronic rnas. eighty-nine percent of the known loci are found to be transcribed in prokaryotic cells prior or subsequent to infection of the host, while 75% of all protein-coding loci are represented in the polyadenylated transcriptomes of human host cells.  dual 3’seq was alternatively coupled to mace to assess the advantages and drawbacks of a library preparation procedure that allows for sequencing of longer fragments. additionally, the identified expression patterns of both organisms were validated by qrt-pcr using three independent biological replicates, which confirmed that relb along with nfkb and nfkb are involved in the initial immune response of epithelial cells after infection with s. enterica typhimurium. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords dual 3’seqdeepsupersagemacetag-basedsimultaneousgenome-widegene expression profilinghost-pathogen interactiontranscriptomesalmonella enterica typhimurium strain sl1344issue-copyright-statement© the author 2015 