if (!requireNamespace("BiocManager", quietly = TRUE))

  install.packages("BiocManager")

BiocManager::install("DESeq2")

#load the DESeq2 package

library(DESeq2)

#load RNA-seq count data

counts <- read.delim("GSE103584_R01_NSCLC_RNAseq.txt", header=TRUE, row.names=1)

#replace NAs with zeros

counts[is.na(counts)] <- 0

#apply the integer conversion to all columns except the first

counts[, -1] <- lapply(counts[, -1], as.integer)

counts$`R01.023` <- as.integer(counts$`R01.023`)

counts$`R01.024` <- as.integer(counts$`R01.024`)

#read in metadata (sample information)

colData <- read.delim("sample_metadata.txt", header=TRUE, row.names=1)

rownames(colData) <- gsub("-", ".", rownames(colData))

#then we use these row names as the column names for counts

colnames(counts) <- rownames(colData)

#Convert the treatment column in colData to a factor

colData$treatment <- as.factor(colData$treatment)

#check for NA values in the counts matrix

if (any(is.na(counts))) {

  cat("NA values found in the counts matrix.\n")

  # Replace NA values with zeros

  counts[is.na(counts)] <- 0

}

#confirm that there are no more NA values

if (any(is.na(counts))) {

  cat("NA values are still present after replacement.\n")

} else {

  cat("No NA values are present. Proceeding with DESeq2.\n")

}

#proceed with DESeq2 if no NAs are present

dds <- DESeqDataSetFromMatrix(countData = counts, colData = colData, design = ~ treatment)

#run the DESeq pipeline

dds <- DESeq(dds)

#running the differential analysis

res <- results(dds, contrast=c("treatment", "male", "female"))

#ordering results by significance

resOrdered <- res[order(res$pvalue),]

#exporting results to csv

write.csv(as.data.frame(resOrdered), file="DESeq2_results.csv", row.names=TRUE)