Turn Data into Meaningful Information

Make it work, then make it fast Atlantbh had the opportunity to implement a back-end software solution for an insurance company. The client required a data processing engine able to transform different sizes of data varying from a few kilobytes to dozens of gigabytes, using the same mechanism. Having a unified engine which would have the same transformation logic both for the files containing 2 records (small data sets) and 20 million records (big data sets) was important to achieve simplicity. From the functional point of view, following requirements had to be addressed: To start a processing flow, user needs to provide the input data, along with the configuration data set. Both are provided using the predefined schema and are validated by the system. Once triggered, each processing flow contains multiple logical steps, each of them transforming the data. Processing flow is dynamic and driven by the configuration. Complexity is very different for each step (from simple operations to application of very extensive business rules). Files given as input and produced as output in each step, often need to be in a different format and sometimes the number of files cannot be predicted in advance. There was also set of non-functional requirements: Performance requirements are given for both small and big data sets and defined as highly important aspect for business. For example, small data sets such as files up to 10k records need to be processed in a range of minutes and should never be blocked by the big ones. There is also a minimum number of processes that should run in parallel, without affecting performance, etc. Different customers should be able to use a system and have full control over their data. They shouldn’t be affected by other users by any means, and access to their data should be secured. It should be possible to control priorities for the different processing flow types. Faced with the combination of business and technical constraints, we needed to come up with a right architecture that will be both easy to extend and maintain. On top of this, a working solution was needed within a short timeframe, in order to enable testing and verification of complex business rules. Challenges Challenge 1: Technology which will support required functionalities needs to be chosen Solution: The choice of technology to be used for data processing and transformation was already made by the client. They decided on Apache Spark, considering that it can support the scalability required by the system. This way the development team is limited by the technology they can use, but are free to use it in a way that will support all business cases – known at the time and possibly needed in the future. Challenge 2: Define how Apache Spark will be used to support all business cases. Possibilities and Solution: To make an adequate decision on how to design the architecture and implement certain functionalities using Spark in a way that would support all business cases, two approaches were considered: usage of Thrift Server, and writing custom Spark jobs using Scala, submitted for processing via custom Process Manager. 1) The usage of Thrift server seamed simpler for implementation but had limitations that would complicate things in the long run With Thrift Server, the main focus would be primarily on describing how the actual data transformations would be performed, as all the responsibilities related to the actual technical usage of the Spark cluster would have been shifted to the Thrift Server. One processing context is used and reused. Jobs would be performed over the same context in case of a small data set. Ensuring that any unnecessary overhead for preparation of processing was removed. On the other hand, the thrift server available inside the used platform seriously limited the flexibility. Only SQL language for transformations was allowed – seriously limiting the transformation capabilities for more complex business cases.  Data IO was limited to the same data format and no option of fine-grained control over resources or parallelism configuration, which is integral for fine-tuning the transformation processing for different sizes of data. 2) Writing custom Spark jobs using Scala was a more complex solution The benefit was that it gave the development team the needed flexibility and freedom in designing the jobs. It required more expertise on the low level of abstraction and would shift the responsibility of resource allocation for specific jobs on the development team. But this way, each processing job would be independent and more flexibility could be achieved, since Scala does not have as many limits as SQL and it can be combined with other languages this way. Considering both options but having in mind the requirements given in the introduction, it was decided that we would go with the second solution. More flexible solution was chosen to ensure that new problems and constraints will not emerge along the way and the team had the expertise to face the challenge. Challenge 3: Meet performance requirements despite the different challenges Moving further, the chosen solution showed certain fallbacks related to the performance. Spark was running on the Hadoop PaaS platform, which used YARN as a resource manager. It meant that each Spark job required the following steps to be executed (List 1*): job acceptance, resource allocation, spinning up of a new context, read resource(s) from HDFS, execution, write output resource(s) to HDFS and shut down used context. This brought unnecessary overhead to each job, which was prominently visible in case of small data sets, where supporting steps could have taken more time than the main execution step. Solution: Due to challenges faced and time constraints given, we agreed with the client to adopt the approach: make it work, then make it fast. We enabled the client to perform business rules validation, while the team worked on the optimizations to make performances better. Optimizations were a multipart process, but we would emphasize three parts that made the most impact. 1) Overhead reduced by chaining smaller jobs into large one Initially, each step was implemented separately, emphasizing the issue with the YARN process manager and stage overheads – going through the same stages (List 1*) for each part of the transformation process. Therefore, we decided to chain all jobs (executions) into one large job, where input and output operations were performed just once at the beginning and the end of the encompassing job instead of each individual step. This required complex framework outside of the jobs to configure, correctly orchestrate and take the burden of reading/writing operations instead of Spark jobs. 2) “Fast-lane” introduced for the small data sets On a larger scale, outside of a single data transformation flow, there were performance overheads with running small data sets on a full-blown actual Spark cluster. The main issue was that, the overheads coming from Spark resource scheduling, contact initialization and driver – executors communication were not justifiable for relatively small job executions which came with smaller datasets. Change has been made to run small data sets in existing pre-configured local Spark contexts, without changing or separating implementation of the business rules. With the introduction of a “fast-lane” for running smaller spark jobs on an already pre-existing context, the overhead from the usual spark application spin-up has been removed – effectively reducing the processing time of small jobs to only the time taken to actually execute the business logic. Large datasets continued to be executed on the Spark cluster, while small jobs would either go to the Spark Local Context or the cluster, depending on client requirements. 3) Maximized throughput with resource utilization configurations For scenarios where both large and small jobs are running on the cluster, there is the natural problem of properly utilizing the available resources to maximize the average throughput between job sizes. Therefore, each job’s resource utilization was dynamically configured based on the dataset size that needed to be processed. A queuing mechanism has been established, with proper job size and resource consumption configurations defined, based on the underlying cluster. Additional factor which we needed to take care of is to allow parallel executions of small and large jobs, but making sure that large jobs would not block the cluster by using all available resources. In the end, these Spark optimizations played a large role in achieving our final goal, meeting customer performance requirements! Conclusion When faced with a short deadline, it may seem more desirable to use technologies, tools or frameworks which require no or limited deep-dive knowledge about how they work under the surface. Nevertheless, we were aware that in this case it wouldn’t be the right approach. We considered not just the actual requirements which were sufficiently complex, but the possibility of them being significantly extended in the future which is very common for the type of the system we worked on. Therefore, we have decided to go for something more complex and flexible, which in the short run may have seemed like a bad decision, but in the long run paid large dividends. Two important factors influenced choosing this path. First, we had the personnel and expertise to be able to quickly adopt all the necessary knowledge. Secondly, we were able to identify with our client what needs to be done to make it ‘just work on time’ and then focused on the ‘make it fast’ part. This case study was written by Atlanters: Nejra Suljić, Irma Smajić, Dino Prašo, and Ridvan Appa Bugis. If you enjoyed it, share it with your friends and colleagues – or read more blogs about Big Data!

A story about adding a new column to a table with 600 million records The Challenge Picture this. You have a task to implement data loading from a large PostgreSQL database table to Snowflake through a data pipeline. Simple enough, right? Just do it incrementally, using the appropriate column as a bookmark to track how much data has already been processed. The problem is that the table has 600 million records in production and lacks a column suitable for the bookmark. This was exactly the situation we faced recently. We needed to: Add an updated_at column to track changes. Backfill existing records with updated_at = created_at Add an index on the updated_at column for efficient querying (essential for bookmark-based incremental loads) Do all of this without affecting production traffic. Why Traditional Approaches Don’t Work Even though created_at was there, it didn’t fit. The updated_at timestamp was the most suitable choice for the bookmark because a single record can be updated over time, and we need to capture those updates.  The naive approach would be to create a new migration with an ALTER TABLE statement, followed by an UPDATE to populate the column, and then CREATE INDEX: ALTER TABLE tableADD COLUMN updated_at TIMESTAMP DEFAULT now();UPDATE tableSET updated_at = created_at;CREATE INDEX idx_table_updated_atON table(updated_at); The whole migration would be wrapped inside a transaction by default.  But what happens next? Adding a column with a default value updates each existing row. This can take hours and locks all rows until it’s done. An UPDATE without batching creates a large transaction that locks all rows for an extended period. Even without the transaction, the lock stays until the statement finishes. This can take a while. So, moving the statement to a separate migration won’t help much. A CREATE INDEX operation gets an exclusive lock on the whole table. This stops all writes and reads for hours.  Any of these operations could bring production to its knees. The table size is the real issue. With less data, all these operations would finish quickly and not disrupt the system. Locking becomes a problem when the operation takes too long. The Zero-Downtime Strategy Here’s how we addressed the issue using PostgreSQL’s features and strategic batching. Step 1: Add the Column without a default value We began by adding the column without immediately populating it. In PostgreSQL, adding a nullable column without a default value is almost instantaneous: ALTER TABLE table ADD COLUMN updated_at TIMESTAMP; Step 2: Batch Update Existing Records Rather than updating all 600 million records in a single transaction, we processed them in smaller chunks of 10,000 records. This method: Keeps transactions small and efficient Allows other operations to continue between batches Prevents lock escalation DO ‘DECLARE chunk_size INT := 10000; start_id INT; max_id INT;BEGIN SELECT MIN(id), MAX(id) INTO start_id, max_id FROM table; WHILE start_id <= max_id LOOP EXECUTE format( ”UPDATE table SET updated_at = created_at WHERE id >= %s AND id < %s”, start_id, start_id + chunk_size ); start_id := start_id + chunk_size; END LOOP;END’; Step 3: Create Index Concurrently PostgreSQL’s CREATE INDEX CONCURRENTLY is a valuable option in these scenarios.  While index creation is slower than usual, it offers greater safety for production environments. Unlike regular index creation, it: Does not block reads and writes Performs multiple table scans to build the index Takes longer but ensures full availability CREATE INDEX CONCURRENTLY idx_table_updated_at ON table (updated_at); Important note: This command must be executed outside of a transaction block. Running it within a transaction would negate its benefits and cause blocking. The complete migration script now looks like this: ALTER TABLE table ADD COLUMN updated_at TIMESTAMP;DO ‘DECLAREchunk_size INT := 10000;start_id INT;max_id INT;BEGINSELECT MIN(id), MAX(id) INTO start_id, max_id FROM table;WHILE start_id <= max_id LOOPEXECUTE format(”UPDATE table SET updated_at = created_at WHERE id >= %s AND id < %s”,start_id, start_id + chunk_size);start_id := start_id + chunk_size;END LOOP;END’;CREATE INDEX CONCURRENTLY idx_table_updated_at ON table (updated_at); Implementation with Liquibase As our backend is built on Spring Boot, we utilize Liquibase for database migrations. Below is a complete migration file that reflects what we implemented in the production environment: <?xml version=”1.1″ encoding=”UTF-8″ standalone=”no”?><databaseChangeLog xmlns=”http://www.liquibase.org/xml/ns/dbchangelog” xmlns:xsi=”http://www.w3.org/2001/XMLSchema-instance” xsi:schemaLocation=”http://www.liquibase.org/xml/ns/dbchangelog http://www.liquibase.org/xml/ns/dbchangelog/dbchangelog-3.5.xsd”> <changeSet id=”add-table-updated_at” author=”amahovac” context=”post-release”> <addColumn tableName=”table”> <column name=”updated_at” type=”TIMESTAMP WITH TIME ZONE”/> </addColumn> </changeSet> <changeSet id=”update-table-updated_at” author=”amahovac” runInTransaction=”false” context=”post-release”> <!– Batch update logic here –> </changeSet> <changeSet id=”create-concurrent-index-table-updated_at” runInTransaction=”false” author=”amahovac” context=”post-release”> <!– Concurrent index creation here –> </changeSet></databaseChangeLog> Key Details Separate migrations We organized operations into distinct migrations (changesets) to improve control over each. Running migrations outside the transaction This is crucial for concurrent operations. If not managed, Liquibase would encompass everything in a transaction, leading to the locks we aim to avoid. Adding runInTransaction=”false” prevents this issue. Skipping migrations during deployment We used context=”post-release” and adjusted the Spring Boot configuration files to set spring.liquibase.contexts=!post-release before deployment.  This way, migrations won’t run during deployment and potentially time out our CodeDeploy process. Running it manually After deployment, we executed these migrations manually during a low-traffic period to monitor performance impacts and minimize risk. The Result The entire operation was completed successfully with zero downtime. The production application continued serving requests throughout the 5+ hour process. Our data pipeline now efficiently syncs the table incrementally to Snowflake.  Lessons Learned Test with production-sized data: Our initial approach performed well in the development environment with 1 million rows, but required significant adjustments for the production database. Always consider data volume and test with realistic datasets. Monitor during execution: We kept an eye on CPU, I/O, and lock metrics throughout the process. Document everything: Future team members need to understand the nuances of these migrations. Our inline comments saved confusion later. Always include created_at and updated_at columns when creating a table: You never know when they might be needed. If we had included them from the start, we could have avoided all these complications. Avoid unnecessary locking: Consider concurrent solutions and batch large updates to minimize lock escalation. Run DDL operations with caution: If possible, execute long-running migrations outside of transactions and separate them from deployment.  Working with large production database tables requires patience, careful planning, and an awareness of the impact of each operation. The extra effort to ensure zero downtime is always worthwhile when users don’t even notice you just modified 600 million records.

RNA sequencing (RNA-seq) has become an indispensable tool for studying transcriptomes. By quantifying RNA molecules across conditions, it enables researchers to explore gene expression patterns, discover novel transcripts, and investigate molecular mechanisms underlying various biological phenomena, including diseases. RNA-seq starts with extracting RNA from cells, sequencing it, and producing millions of short reads. These reads are then aligned to a reference genome to produce a matrix of counts representing the number of reads per gene. However, raw counts alone are not sufficient for meaningful biological insights due to technical and biological variations. This is where Differential Expression Analysis (DEA) comes into play, quantifying and statistically validating differences in gene expression between experimental groups. There are many tools and bioinformatics programs used for the RNA-seq differential expression analysis. Why Use Bioconductor for RNA-Seq Analysis? Bioconductor is an R-based ecosystem specifically designed for bioinformatics. It includes packages like DESeq2, which handle the complexities of RNA-seq analysis, such as normalization, model fitting, and visualization. Bioconductor emphasizes reproducibility, a cornerstone of good scientific practice, while providing flexibility for diverse experimental designs. The statistical rigor of its tools ensures that results are both accurate and interpretable, making it an ideal choice for RNA-seq workflows. Introduction to the RNA-Seq Workflow Before diving into the differential expression analysis, it is essential to understand the broader RNA-Seq workflow and its foundational steps. RNA-Seq is not just about analyzing gene expression differences; it begins with careful experimental design and follows a structured pipeline to ensure reliable results. A robust experimental design is crucial, as errors at this stage can compromise the entire analysis. Biological replicates, typically three or more, are essential for capturing natural variability and improving the robustness of statistical tests. Minimizing batch effects and avoiding confounding variables, such as sex or age differences, is equally important to ensure that observed changes in gene expression reflect true biological differences rather than experimental artifacts. The workflow proceeds with sample preparation, where RNA is isolated, and non-essential components like ribosomal RNA are removed. The RNA is then reverse-transcribed into cDNA, fragmented, and prepared into libraries for sequencing. Sequencing generates millions of short reads, which are output as FASTQ files, representing the raw data for analysis. Quality control is performed on these files to assess sequencing accuracy and check for contamination or technical issues. Next, the reads are aligned to a reference genome, mapping them to their most likely origin in the genome. This alignment is particularly complex in RNA-Seq, as mRNA consists only of exons, and many reads span intronic boundaries, necessitating tools capable of aligning across introns. Once aligned, reads are quantified to produce a count matrix, a table where rows represent genes and columns represent samples. This count matrix forms the input for statistical analyses, where the primary goal is to determine whether gene expression differences between conditions are statistically significant, accounting for natural variability within sample groups. The results of the differential expression analysis include log2 fold changes and adjusted p-values for each gene, identifying those with significant expression differences. This process not only pinpoints individual genes of interest but also highlights broader patterns and pathways that may drive the observed biological phenomena. Bioconductor workflow Step 1: Data Preparation and Preprocessing Let’s try to perform DEA with Bioconductor using a publicly available RNA-seq dataset for fibrosis in Mus musculus (house mouse)*. Fibrosis is a pathological condition in which connective tissue replaces normal parenchymal tissue with connective tissue during wound healing, which ultimately reduces organ function and can cause failure and death**. The dataset used for this analysis contains raw RNA counts for each gene in 4 samples with induced fibrosis condition in mice through overexpression of the SMOC2 gene, which has been found to play one of the key roles in fibrosis condition, and 3 normal samples without the pathology.  The first step would be installing all of the packages needed for this analysis and then loading our raw count matrix. {r} # Install BiocManager (for Bioconductor packages like DESeq2) install.packages(“BiocManager”) # Install DESeq2 BiocManager::install(“DESeq2”) A raw count matrix from RNA-Seq is a table that contains the raw, unnormalized counts of sequencing reads that map to genomic features (e.g., genes, transcripts, or exons). These counts represent the initial data obtained after read alignment and quantification. After loading our raw data, it should look something like this: Step 2: Metadata Creation and Alignment Metadata describes the experimental conditions for each sample. In this case, metadata includes information about treatment groups “fibrosis” and “normal.” Proper metadata creation ensures each sample is correctly linked to its condition, forming the basis for statistical comparisons. {r} metadata <- data.frame( treatment = c(“smoc2”, “smoc2”, “smoc2”, “smoc2”, “smoc2”, “smoc2”, “smoc2”), group =c(“fibrosis”, “fibrosis”, “fibrosis”, “fibrosis”, “normal”, “normal”, “normal”) ) rownames(metadata) <- colnames(data) The resulting metadata should look something like this: The alignment between the metadata and the count matrix is verified to prevent mismatches, which could lead to invalid results. This step reflects careful experimental planning, as the metadata directly determines how the statistical model interprets the data. Step 3: Creating the DESeq2 Dataset Using the DESeq2 package, we encapsulate the count data, metadata, and experimental design into a DESeq2 dataset object. The design formula specifies the condition of interest, such as the “group” variable in this case, and accounts for sources of variation. Output from this step can be observed below. {r} dds <- DESeqDataSetFromMatrix( countData = data, colData = metadata, design = ~ group ) This step sets the stage for differential expression analysis by modeling RNA-seq data using the negative binomial distribution. This distribution accommodates the variability typical of RNA-seq experiments, which often exceeds what simpler distributions like the Poisson can handle. Step 4: Data Normalization Raw RNA-seq counts are affected by sequencing depth and library size differences. To address this, DESeq2 calculates size factors for each sample, normalizing the data to make gene expression levels comparable across samples. {r} dds <- estimateSizeFactors(dds) normalized_counts <- counts(dds, normalized = TRUE) Normalization adjusts for biases and ensures that observed differences in gene expression are due to biological, rather than technical, factors. Skipping normalization could lead to inaccurate identification of differentially expressed genes, as systematic biases would confound the analysis. Step 5: Quality Assessment with VST, Heatmaps, and PCA After normalization, it is crucial to assess the data quality and explore relationships between samples. This ensures that biological replicates cluster together and conditions separate as expected. Variance Stabilizing Transformation (VST) RNA-seq data often exhibits a dependence of variance on mean expression, which can obscure patterns in the data. To address this, DESeq2 applies a Variance Stabilizing Transformation (VST) that makes variance approximately constant across expression levels. This transformation enhances the interpretability of downstream clustering and visualizations. {r} vsd <- vst(dds, blind = TRUE) The blind = TRUE argument ensures the transformation is not biased by experimental design, making it ideal for quality assessment. Heatmaps for Sample Clustering Heatmaps provide a visual representation of sample similarity based on pairwise correlation of gene expression. Samples are hierarchically clustered, revealing whether replicates group together and conditions separate clearly. A correlation heatmap shows the pairwise correlation coefficients between samples based on their VST-transformed expression data. These coefficients range from -1 (completely dissimilar) to 1 (perfectly similar). This heatmap reveals how closely related each sample is to the others in terms of overall gene expression. {r} pheatmap(cor(assay(vsd)), main = “Sample Correlation Heatmap”) If replicates cluster tightly and conditions are well-separated, like it is the case in our heatmap above, it indicates good experimental reproducibility. Unexpected clustering patterns or outliers might require further investigation. Distance Heatmap A distance heatmap calculates the pairwise Euclidean distances between samples based on their VST-transformed data. It provides an alternative perspective on sample relationships, focusing on absolute differences in expression profiles rather than relative similarity. {r} pheatmap(as.matrix(dist(t(assay(vsd)))), main = “Sample Distance Heatmap”) In the distance heatmap, smaller distances (dark red) indicate greater similarity between samples. Biological replicates should form tight clusters with small distances, while samples from different conditions should exhibit larger distances (dark blue). This visualization is particularly sensitive to outliers, as samples with unusual expression profiles will appear far from all others. Principal Component Analysis (PCA) PCA is another powerful technique to examine sample similarity. It reduces data dimensionality, with the first principal component (PC1) capturing the largest source of variation. In RNA-seq, PCA often reveals whether the experimental condition drives variation across samples. {r} plotPCA(vsd, intgroup = “group”) In the PCA plot, biological replicates should cluster, while different conditions should separate along one of the principal components. Any failure to see separation might indicate confounding factors or subtle biological differences. We can see that in our sample groups, normal and fibrosis, separate well on PC1. This means that our condition corresponds to PC1, which represents 93% of the variance in our data, while 2% is explained by PC2. This is great since it seems that a lot of the variation in gene expression in the dataset can likely be explained by the differences between sample groups.  Step 6: Differential Expression Analysis The differential expression analysis aims to identify genes with significant expression differences between conditions, such as fibrosis and normal. DESeq2 fits the count data to a negative binomial model and performs statistical tests to compute log2 fold changes and adjusted p-values for each gene using the DESeq2()function in Bioconductor. {r} dds <- DESeq(dds) results <- results(dds, contrast = c(“group”, “fibrosis”, “normal”), alpha = 0.05) The alpha parameter represents the significance threshold for determining which genes are considered “differentially expressed” (DE). An alpha of 0.05 means we are allowing up to a 5% chance of identifying a gene as DE when it is actually not DE (a false positive). This threshold corresponds to a p-value cutoff of 0.05. Genes with adjusted p-values (padj) below 0.05 are considered statistically significant. Exploring Gene Variability with Dispersion Plot The dispersion plot is an essential diagnostic tool for RNA-seq analysis. It shows how variability (dispersion) in gene expression relates to mean expression levels. {r} plotDispEsts(dds) Black Dots: Raw dispersion estimates for each gene. Red Line: Expected dispersion trend based on the mean expression. Blue Dots: Shrunken estimates used in the model for better accuracy. A good fit is indicated by points clustering around the red line. Large deviations might signal poor data quality or unaccounted variability. Shrinking dispersion values ensures more reliable differential expression results, especially for genes with low counts. MA Plot Another way to visually provide insight into gene behavior and ensure the reliability of differential expression analysis is an MA plot. The MA plot visualizes the relationship between the mean expression and the log2 fold change for all genes, helping to identify differentially expressed genes. {r} plotMA(results, ylim = c(-8, 8)) Genes with significant expression changes are highlighted, but genes with low counts may show large, unreliable fold changes. In this case, we can apply the shrinkage method. Shrinking reduces exaggerated fold changes for genes with low counts, improving reliability without affecting statistical significance. {r} results_shrunken <- lfcShrink( dds, contrast = c(“group”, “fibrosis”, “normal”), type = “normal”) plotMA(results_shrunken, ylim = c(-8, 8)) The difference between normal and shrunken MA plots can be observed below. Refining Results After identifying differentially expressed genes (DEGs) with DESeq2, we apply thresholds to focus on biologically significant findings to ensure selected genes show both meaningful expression changes and statistical significance. The lfcThreshold specifies the minimum log2 fold change (LFC) that is considered biologically meaningful. A value of 0.32 corresponds to an absolute fold change of approximately 1.25. The alpha parameter sets the significance level for the adjusted p-value (Benjamini-Hochberg False Discovery Rate), and a value of 0.05 means we accept a 5% false discovery rate. These thresholds can be adjusted based on the study’s requirements. {r} results_lfc <- results( dds, contrast = c(“group”, “fibrosis”, “normal”), lfcThreshold = 0.32, alpha = 0.05 ) To add biological context, Ensembl gene IDs are annotated using a very convenient Bioconductor package: the biomaRt: {r} mart <- useMart(“ensembl”, dataset = “mmusculus_gene_ensembl”) annotations <- getBM( attributes =c(“ensembl_gene_id”, “external_gene_name”, “description”), filters = “ensembl_gene_id”, values =significant_results_lfc$ensgene, mart = mart ) Our significant differentially expressed genes with their biological descriptions and roles can be clearly observed now: Step 7: Visualization of Results Visualizations like volcano plots and heatmaps make it easier to interpret differential expression results. The expression heatmap visualizes the normalized counts of significant genes (in this case, the top 20 were chosen), highlighting patterns of gene expression across samples. {r} pheatmap(sig_norm_counts, color = heat_colors, scale =”row”) Genes and samples with similar expression profiles group together, and it can be seen how significant genes are expressed differently between sample groups. A Z-score standardizes the gene expression values for each gene (row) to highlight relative expression patterns. Negative values (e.g., -2) show that expression is lower than average for that gene. Positive values (e.g., 2) show that expression is higher than average for that gene. Another way to visualize the expression levels of the most significant genes across all samples is the expression plot. {r} ggplot(top_20, aes(x = ensgene, y = normalized_counts, color = group)) + geom_point(size = 3, alpha = 0.8) + scale_y_log10() + labs( title = “Expression of Top 20 Genes”, x = “Gene ID”, y = “Normalized Counts (log10 scale)”, color = “Condition” ) + theme_minimal() + theme( axis.text.x = element_text(angle = 45, hjust = 1), plot.title =element_text(hjust = 0.5) ) The plot highlights the variation in the expression of these genes between the sample groups, and it is easy to observe upregulated or downregulated genes between the conditions. Conclusion and Next Steps These steps show the basic workflow for RNA-seq differential expression analysis using DESeq2, from preprocessing data to identifying significant differentially expressed genes (DEGs). While these results are valuable, they are just the beginning. The next steps involve exploring and validating the findings: Research Genes of Interest: Focusing on the most significant DEGs. Investigating their known functions, roles in relevant pathways, or associations with specific conditions using databases like Ensembl or Gene Ontology. Experimental Validation: Confirming the computational findings through laboratory experiments such as qPCR or western blotting. Functional Analysis: Performing pathway enrichment or gene ontology analysis to uncover biological processes and pathways associated with the significant genes. Bioconductor provides tools like clusterProfiler for this. Hypothesis Generation: Using results to form new hypotheses about the biological mechanisms underlying the condition of interest. By following these steps, RNA-seq data can be turned into meaningful biological insights, driving further research and discovery. References *Oakley, F., Gee, L. M., Sheerin, N. S., & Borthwick, L. A. (2019). Implementation of pre-clinical methodologies to study fibrosis and test anti-fibrotic therapy. Current Opinion in Pharmacology, 49, 95–101. https://doi.org/10.1016/j.coph.2019.10.004  **Oakley, F., Gee, L. M., Sheerin, N. S., & Borthwick, L. A. (2019). Implementation of pre-clinical methodologies to study fibrosis and test anti-fibrotic therapy. Current Opinion in Pharmacology, 49, 95–101. https://doi.org/10.1016/j.coph.2019.10.004 Work was inspired by Datacamp’s RNA-Seq with Bioconductor in R course. “Using Bioconductor for RNA-seq Differential Expression Analysis” Tech Bite was brought to you by Zehra Alibašić, Junior Data Analyst at Atlantbh. Tech Bites are tips, tricks, snippets or explanations about various programming technologies and paradigms, which can help engineers with their everyday job.