1 Overview

batch condition media treatment induced exp_11
SK061 1 LB-DMSO LB DMSO no control
SK062 2 LB-DMSO LB DMSO no control
SK063 3 LB-DMSO LB DMSO no control
SK082 1 SCFM-DMSO SCFM DMSO no treatment
SK083 2 SCFM-DMSO SCFM DMSO no treatment
SK084 3 SCFM-DMSO SCFM DMSO no treatment

2 Project variables

##       batch condition media treatment induced    exp_11
## SK061     1   LB-DMSO    LB      DMSO      no   control
## SK062     2   LB-DMSO    LB      DMSO      no   control
## SK063     3   LB-DMSO    LB      DMSO      no   control
## SK082     1 SCFM-DMSO  SCFM      DMSO      no treatment
## SK083     2 SCFM-DMSO  SCFM      DMSO      no treatment
## SK084     3 SCFM-DMSO  SCFM      DMSO      no treatment
## [1] "The control is: LB-DMSO"
## [1] "The treatment is: SCFM-DMSO"
# Load raw count data
rawcounts <- read.csv("../../raw-data/htseq_SK_unfiltered_reverse_strand_count.csv", 
                      row.names = 1, header = TRUE)

# Append a prefix to the batch and run columns to turn to character class
metadata[["batch"]] <- paste0("b", metadata[["batch"]])
metadata[["run"]] <- paste0("r", metadata[["run"]])

# Subset the rawcounts dataframe to only include samples that are present in experiment
samples <- rownames(metadata)
rawcounts <- select(rawcounts, all_of(samples))

# Ensure that the order of samples in metadata matches the order of columns in rawcounts
all(rownames(metadata) == colnames(rawcounts))  # Should return TRUE if order is correct
## [1] TRUE
# Define the design formula for DESeq2, including main condition and batch effect
design_formula <- ~ exp_11 + batch

# Create DESeq2 dataset object from the filtered count matrix and metadata
dds <- DESeqDataSetFromMatrix(
  countData = rawcounts,
  colData = metadata,
  design = design_formula
)
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors

3 Library size per sample

# Sum up the counts for each sample
rawcounts_sums <- colSums(rawcounts)

# Turn it into a data frame
data <- data.frame(
  Category = names(rawcounts_sums),
  Count = as.numeric(rawcounts_sums)
)

# Make sure sample order is consistent
data <- data %>%
  mutate(Category = factor(Category, levels = sort(unique(Category))))

# Plot bar chart of library sizes
p <- ggplot(data, aes(x = Category, y = Count, fill = Category)) +
  geom_bar(stat = "identity") +
  # coord_flip() # flip if too many samples to read labels easily
  labs(y = "Library size (M)", x = NULL) +  # remove x label, rename y axis
  # scale_fill_brewer(palette = "Dark2") +  # optional color palette
  scale_y_continuous(
    labels = scales::label_number(scale = 1e-6),  # show counts in millions
    limits = c(0, (max(data$Count) + 1e6))        # add a bit of space on top
  ) +
  theme(
    panel.grid = element_blank(),
    panel.background = element_blank(),
    axis.text = element_text(size = text.size, angle = 45, hjust = 1),  # slant labels
    axis.title = element_text(size = axis.text.size),
    legend.text = element_text(size = text.size),
    legend.title = element_text(size = axis.text.size),
    panel.border = element_rect(fill = NA),
    legend.position = "none",  # no legend needed since each bar is a different sample
    axis.ticks = element_line(linewidth = 0.2),
    legend.key = element_blank()
  )

# Show the plot
p

# Save the plot as a PDF with today’s date and project name
ggsave(paste0(Sys.Date(), "_library_size_", exp, "_", treatment,".vs.", control, ".pdf"), device = "pdf")
## Saving 7 x 5 in image

4 Expression distribution

# Stack rawcounts into single vector 
dat <- stack(rawcounts)

# Format expression distribution plot with the expression level and density 
expression_dist <- ggplot(dat, aes(x = log2(values+1),  fill=ind, label = ind )) + 
  geom_density(alpha = 0.5, position = "identity") + 
  ylab("Density") + 
  geom_textdensity(hjust = "ymax", vjust = 0,
           text_only = TRUE, text_smoothing = 20) + 
    theme(panel.grid = element_blank(),
        panel.background = element_blank(),
        axis.text = element_text(size = text.size),
        axis.title = element_text(size = axis.text.size),
        legend.text = element_text(size = text.size),
        legend.title = element_text(size = axis.text.size),
        panel.border = element_rect(fill = NA),
        legend.position = "none",
        axis.ticks = element_line(linewidth = 0.2),
        legend.key = element_blank()) +
  coord_cartesian(xlim = c(0, 15))

expression_dist

#save figure
cairo_pdf(paste0(Sys.Date(), "_expression_dist_", exp, "_", treatment, ".vs.", control, ".pdf"))
expression_dist
dev.off() 
## quartz_off_screen 
##                 2

5 Non-zero genes observed

# Prepare data by getting library sizes and non-zero counts
x <- colSums(rawcounts)
y <- colSums(rawcounts != 0) 
df <- data.frame(x,y)
df$ind <- rownames(df)

# Calculate axis limits
x_min <- min(df$x)
x_max <- max(df$x)
y_min <- min(df$y)
y_max <- max(df$y)

# Non-zero plot
nonzero <- ggplot(data = df, aes(x = x, y = y, color = ind, label = ind)) +
  geom_point(size = 3) +
  geom_text(hjust = -0.1, vjust = -0.1, color = "black", size = 5) +
  scale_y_continuous(limits=c(y_min, y_max)) + 
  scale_x_continuous(
    limits = c(x_min-1e3, x_max+9e5),
    labels = scales::comma_format(scale = 1e-6,  accuracy = 1)  # Format labels without decimals
  ) +
  xlab("Millions of reads") + ylab("Number of non-zero genes observed") +
  labs(color = "Sample") + 
    theme(panel.grid = element_blank(),
        panel.background = element_blank(),
        axis.text = element_text(size = text.size),
        axis.title = element_text(size = axis.text.size),
        legend.text = element_text(size = text.size),
        legend.title = element_text(size = axis.text.size),
        panel.border = element_rect(fill = NA),
        legend.position = "none",
        axis.ticks = element_line(linewidth = 0.2),
        legend.key = element_blank()) +
  coord_cartesian()

# Show nonzero plot
nonzero

# Save plot

cairo_pdf(paste0(Sys.Date(), "_non-zero_", exp, "_", treatment, ".vs.", control, ".pdf"))
nonzero
dev.off() 
## quartz_off_screen 
##                 2

6 RPKM and CPM transformations

# Load GFF annotation 
gff_annot <- load_gff_annotations("../../raw-data/paeruginosa_pa14.gff", 
                                  id_col = "gene_id", type = "gene")
## Returning a df with 16 columns and 5979 rows.
rownames(gff_annot) <- gff_annot$gene_id

head(gff_annot)
##               seqnames start  end width strand    source type score phase
## gene1650835 chromosome   483 2027  1545      + PseudoCAP gene    NA     0
## gene1650837 chromosome  2056 3159  1104      + PseudoCAP gene    NA     0
## gene1650839 chromosome  3169 4278  1110      + PseudoCAP gene    NA     0
## gene1650841 chromosome  4275 6695  2421      + PseudoCAP gene    NA     0
## gene1650843 chromosome  7018 7791   774      - PseudoCAP gene    NA     0
## gene1650845 chromosome  7803 8339   537      - PseudoCAP gene    NA     0
##                 gene_id Name         Dbxref      Alias       name Parent locus
## gene1650835 gene1650835      GeneID:4384099 PA14_00010       dnaA         <NA>
## gene1650837 gene1650837      GeneID:4384100 PA14_00020       dnaN         <NA>
## gene1650839 gene1650839      GeneID:4384101 PA14_00030       recF         <NA>
## gene1650841 gene1650841      GeneID:4384102 PA14_00050       gyrB         <NA>
## gene1650843 gene1650843      GeneID:4385186 PA14_00060 PA14_00060         <NA>
## gene1650845 gene1650845      GeneID:4385187 PA14_00070 PA14_00070         <NA>
head(rawcounts)
##             SK061 SK062 SK063 SK082 SK083 SK084
## gene1650835  1555  1680  1858   150   120   137
## gene1650837  1560  1390  1898   338   251   409
## gene1650839   558   436   675   353   264   300
## gene1650841  2492  2352  2588   266   226   196
## gene1650843   285   247   296   135   128   159
## gene1650845   129   106   145   125    70   101
dim(gff_annot)
## [1] 5979   16
dim(rawcounts)
## [1] 5979    6
summary(rownames(gff_annot) == rownames(rawcounts))
##    Mode   FALSE    TRUE 
## logical     932    5047
gff_annot <- gff_annot[rownames(rawcounts),]
summary(rownames(gff_annot) == rownames(rawcounts))
##    Mode    TRUE 
## logical    5979
# Create an expression object with counts, metadata, and gene info
expt <- create_expt(count_dataframe = rawcounts,
                    metadata = metadata, 
                    gene_info = gff_annot)
## Reading the sample metadata.
## Did not find the column: sampleid.
## The rownames do not appear numeric, using them.
## The sample definitions comprises: 6 rows(samples) and 8 columns(metadata fields).
## Matched 5979 annotations and counts.
## Bringing together the count matrix and gene information.
## Some annotations were lost in merging, setting them to 'undefined'.
## Saving the expressionset to 'expt.rda'.
## The final expressionset has 5979 features and 6 samples.
# Normalize counts log2 CPM
cpm_normalized_exp <- normalize_expt(expt, 
                                     transform = "log2", 
                                     norm = "raw", 
                                     convert = "cpm", 
                                     row_min = "10")
## transform_counts: Found 409 values equal to 0, adding 1 to the matrix.
####### Calculate RPKM ########

# Grab gene lengths from GFF (should be in base pairs)
gene_lengths <- gff_annot$width  
names(gene_lengths) <- rownames(gff_annot)  # set gene names

# Make sure gene_lengths are in the same order as rawcounts rows
gene_lengths <- gene_lengths[rownames(rawcounts)]

# Make DGEList object for edgeR
dge <- DGEList(counts = rawcounts)

# Calculate RPKM (reads per kilobase per million)
rpkm <- rpkm(dge, gene.length = gene_lengths)

# Wrap RPKM in an experiment object too
expt_rpkm <- create_expt(count_dataframe = rpkm,
                         metadata = metadata, 
                         gene_info = gff_annot)
## Reading the sample metadata.
## Did not find the column: sampleid.
## The rownames do not appear numeric, using them.
## The sample definitions comprises: 6 rows(samples) and 8 columns(metadata fields).
## Matched 5979 annotations and counts.
## Bringing together the count matrix and gene information.
## Some annotations were lost in merging, setting them to 'undefined'.
## Saving the expressionset to 'expt.rda'.
## The final expressionset has 5979 features and 6 samples.
# log2 transform the RPKM (again, no normalization beyond log)
rpkm_normalized_exp <- normalize_expt(expt_rpkm, 
                                      transform = "log2", 
                                      norm = "raw")
## transform_counts: Found 409 values equal to 0, adding 1 to the matrix.
# Get the actual expression matrices from the experiment objects
rpkm <- exprs(rpkm_normalized_exp)
cpm <- exprs(cpm_normalized_exp)

# Merge CPM and RPKM matrices with gene annotation info
cpm_annotated <- merge(gff_annot, cpm, by = 0)
rpkm_annotated <- merge(gff_annot, rpkm, by = 0)

6.1 Heatmap of alginate genes

This generates a small heatmap of the alginate genes across samples to visualize expression.

# Loop through CPM and RPKM to plot both heatmaps
for (type in c("cpm", "rpkm")) {
  
  # Define gene list
  alg_operon_genes <- c(
    "algA", "algF", "algJ", "algI", "algL", "algX", "algG", "algE",
    "algK", "alg44", "alg8", "algD", "mucC", "mucB", "algU", "algW",
    "algP", "algQ", "algR", "algZ", "algC", "algB"
  )

  # Get the right data per type of transformation
  data_annot <- get(paste0(type, "_annotated"))
  
  # Fix rownames using gene name (make unique)
  rownames(data_annot) <- make.names(data_annot$name, unique = TRUE)
  
  # Drop undesired metadata columns
  cols_to_exclude <- c("Row.names","seqnames","start","end","width",
                       "strand","source","type","score","phase",
                       "gene_id","Name","Dbxref","Alias","name",
                       "Parent","locus")
  keep <- !names(data_annot) %in% cols_to_exclude
  expr_data <- data_annot[, keep]
  
  # Extract matrix for alg genes
  alg_genes <- expr_data[alg_operon_genes, , drop = FALSE]
  mat <- alg_genes - rowMeans(alg_genes)
  
  # Pull sample annotations
  anno <- as.data.frame(metadata[, c("condition", "induced")])
  colnames(anno) <- c("Condition", "Induced")
  
  # Color settings for annotation
 # ann_colors <- list(
 #  Condition = c("LB-DMSO" = "#D95F02", "LB-IPTG" = "#66A61E"),
 #  Induced = c("no" = "lightgrey", "yes" = "gray27")
 # )
  
  # Scale data by row (Z-score per gene)
  scaled_mat <- t(scale(t(mat)))

  # Build heatmap
  ht <- Heatmap(
    scaled_mat,
    name = "Z-score",
    top_annotation = HeatmapAnnotation(df = anno#, 
                                       #col = ann_colors
                                       ),
    cluster_rows = TRUE,
    cluster_columns = TRUE,
    row_names_side = "left",
    column_names_side = "bottom",
    column_names_rot = 45,
    row_names_gp = gpar(fontsize = 10),
    column_names_gp = gpar(fontsize = 10),
    show_row_dend = TRUE,
    show_column_dend = TRUE
  )
  
  # Draw and save
  draw(ht)
  
  cairo_pdf(paste0(Sys.Date(), "_alg_heatmap_", type, "_", exp, "_", treatment, ".vs.", control, ".pdf"))
  draw(ht)
  dev.off()
}

7 DESeq

# Create DESeqDataSet
dds <- DESeqDataSetFromMatrix(
  countData = rawcounts,
  colData = metadata,
  design = design_formula
)
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
# Check number of genes before filtering
print("Number of genes before filtering:")
## [1] "Number of genes before filtering:"
nrow(dds)
## [1] 5979
# Filter: keep genes with counts >=10 in at least 2 samples
print("Filtering genes with counts >=10 in at least 2 samples")
## [1] "Filtering genes with counts >=10 in at least 2 samples"
keep <- rowSums(counts(dds) >= 10) >= 2 
dds_filtered <- dds[keep, ]

# Check number of genes after filtering
print("Number of genes after filtering:")
## [1] "Number of genes after filtering:"
nrow(dds_filtered)
## [1] 5746
# Run DESeq on the filtered dataset
dds_final <- DESeq(dds_filtered)
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
# Show available results
resultsNames(dds_final)
## [1] "Intercept"                   "exp_11_treatment_vs_control"
## [3] "batch_b2_vs_b1"              "batch_b3_vs_b1"
# Plot dispersion estimates
plotDispEsts(dds_final)

8 Save DE tables

res <- results(dds_final, contrast = c("exp_11", "treatment", "control"))
resdf <- as.data.frame(res)
resdf$gene_id <- rownames(resdf)

cpm <- as.data.frame(cpm)
rpkm <- as.data.frame(rpkm)

# append type of value to the dataframe   
colnames(rpkm) <- paste0(colnames(rpkm), "_rpkm")
colnames(cpm) <- paste0(colnames(cpm), "_cpm")
rawcounts2 <- rawcounts
colnames(rawcounts2) <- paste0(colnames(rawcounts2), "_rawcounts")

# ensure column gene_id is there
rpkm$gene_id <- rownames(rpkm)
cpm$gene_id <- rownames(cpm)
rawcounts2$gene_id <- rownames(rawcounts2)

gff_annot <- gff_annot[rownames(resdf),]

dim(resdf)
## [1] 5746    7
# merge with siggenes
resdf <- full_join(resdf, gff_annot, by = "gene_id", copy = FALSE)
head(resdf)
##    baseMean log2FoldChange     lfcSE        stat       pvalue         padj
## 1  786.5229     -3.1219412 0.2140327 -14.5862786 3.434268e-48 6.090526e-47
## 2  865.7989     -1.7706474 0.1993684  -8.8812834 6.609336e-19 4.006038e-18
## 3  411.3546     -0.3275486 0.2205452  -1.4851766 1.374970e-01 1.833081e-01
## 4 1164.6519     -2.9177741 0.1998135 -14.6024908 2.707692e-48 4.831800e-47
## 5  198.0543     -0.4497606 0.2497138  -1.8011040 7.168650e-02 1.021604e-01
## 6  111.0326      0.1408506 0.3099986   0.4543588 6.495706e-01 7.074361e-01
##       gene_id   seqnames start  end width strand    source type score phase
## 1 gene1650835 chromosome   483 2027  1545      + PseudoCAP gene    NA     0
## 2 gene1650837 chromosome  2056 3159  1104      + PseudoCAP gene    NA     0
## 3 gene1650839 chromosome  3169 4278  1110      + PseudoCAP gene    NA     0
## 4 gene1650841 chromosome  4275 6695  2421      + PseudoCAP gene    NA     0
## 5 gene1650843 chromosome  7018 7791   774      - PseudoCAP gene    NA     0
## 6 gene1650845 chromosome  7803 8339   537      - PseudoCAP gene    NA     0
##   Name         Dbxref      Alias       name Parent locus
## 1      GeneID:4384099 PA14_00010       dnaA         <NA>
## 2      GeneID:4384100 PA14_00020       dnaN         <NA>
## 3      GeneID:4384101 PA14_00030       recF         <NA>
## 4      GeneID:4384102 PA14_00050       gyrB         <NA>
## 5      GeneID:4385186 PA14_00060 PA14_00060         <NA>
## 6      GeneID:4385187 PA14_00070 PA14_00070         <NA>
pa_go <- read.csv("../../raw-data/208963_microbes.online.csv", header = TRUE)
head(pa_go)
##   locusId   accession        GI scaffoldId start stop strand      Alias
## 1 2194572 YP_788156.1 116053721       4582   483 2027      + PA14_00010
## 2 2194573 YP_788157.1 116053722       4582  2056 3159      + PA14_00020
## 3 2194574 YP_788158.1 116053723       4582  3169 4278      + PA14_00030
## 4 2194575 YP_788159.1 116053724       4582  4275 6695      + PA14_00050
## 5 2194576 YP_788160.1 116053725       4582  7791 7018      - PA14_00060
## 6 2194577 YP_788161.1 116053726       4582  8339 7803      - PA14_00070
##         name                                                  desc     COG
## 1       dnaA chromosomal replication initiator protein DnaA (NCBI)  COG593
## 2       dnaN                 DNA polymerase III, beta chain (NCBI)  COG592
## 3       recF        DNA replication and repair protein RecF (NCBI) COG1195
## 4       gyrB                           DNA gyrase subunit B (NCBI)  COG187
## 5 PA14_00060                       putative acyltransferase (NCBI)  COG204
## 6 PA14_00070                putative histidinol-phosphatase (NCBI)  COG241
##   COGFun
## 1      L
## 2      L
## 3      L
## 4      L
## 5      I
## 6      E
##                                                                    COGDesc
## 1                            ATPase involved in DNA replication initiation
## 2                      DNA polymerase sliding clamp subunit (PCNA homolog)
## 3                         Recombinational DNA repair ATPase (RecF pathway)
## 4 Type IIA topoisomerase (DNA gyrase/topo II, topoisomerase IV), B subunit
## 5                           1-acyl-sn-glycerol-3-phosphate acyltransferase
## 6                          Histidinol phosphatase and related phosphatases
##                                                                                  TIGRFam
## 1                        TIGR00362 chromosomal replication initiator protein DnaA [dnaA]
## 2                                      TIGR00663 DNA polymerase III, beta subunit [dnaN]
## 3                               TIGR00611 DNA replication and repair protein RecF [recF]
## 4                                                 TIGR01059 DNA gyrase, B subunit [gyrB]
## 5                                                                                       
## 6 TIGR01656 histidinol-phosphate phosphatase domain,TIGR01662 HAD hydrolase, family IIIA
##                                                   TIGRRoles
## 1 DNA metabolism:DNA replication, recombination, and repair
## 2 DNA metabolism:DNA replication, recombination, and repair
## 3 DNA metabolism:DNA replication, recombination, and repair
## 4 DNA metabolism:DNA replication, recombination, and repair
## 5                                                          
## 6           Unknown function:Enzymes of unknown specificity
##                                                                                                                                                                                           GO
## 1                                                                                                                                     GO:0006270,GO:0006275,GO:0003688,GO:0017111,GO:0005524
## 2 GO:0006260,GO:0003677,GO:0003893,GO:0008408,GO:0016449,GO:0019984,GO:0003889,GO:0003894,GO:0015999,GO:0016450,GO:0003890,GO:0003895,GO:0016000,GO:0016451,GO:0003891,GO:0016448,GO:0016452
## 3                                                                                                                                     GO:0006281,GO:0005694,GO:0005524,GO:0017111,GO:0003697
## 4                                                                                                                                     GO:0006304,GO:0006265,GO:0005694,GO:0003918,GO:0005524
## 5                                                                                                                                                                      GO:0008152,GO:0003841
## 6                                                                                                                                                                      GO:0000105,GO:0004401
##         EC                                        ECDesc
## 1                                                       
## 2  2.7.7.7                  DNA-directed DNA polymerase.
## 3                                                       
## 4 5.99.1.3          DNA topoisomerase (ATP-hydrolyzing).
## 5 2.3.1.51 1-acylglycerol-3-phosphate O-acyltransferase.
## 6  3.1.3.-
col_to_remove <- intersect(colnames(pa_go),colnames(resdf))
col_to_remove <- setdiff(col_to_remove, "Alias") # returns all except the second vector

pa_go <- pa_go %>% select(-all_of(col_to_remove))

resdf <- full_join(resdf, rpkm, by = "gene_id", copy = FALSE)
resdf <- full_join(resdf, cpm, by = "gene_id", copy = FALSE)
resdf <- full_join(resdf, rawcounts2, by = "gene_id", copy = FALSE)
#resdf <- full_join(resdf, pa_go, by = "Alias", copy = FALSE)
resdf <- merge(resdf, pa_go, by = "Alias", copy = FALSE)

head(resdf)
##        Alias  baseMean log2FoldChange     lfcSE        stat       pvalue
## 1 PA14_00010  786.5229     -3.1219412 0.2140327 -14.5862786 3.434268e-48
## 2 PA14_00020  865.7989     -1.7706474 0.1993684  -8.8812834 6.609336e-19
## 3 PA14_00030  411.3546     -0.3275486 0.2205452  -1.4851766 1.374970e-01
## 4 PA14_00050 1164.6519     -2.9177741 0.1998135 -14.6024908 2.707692e-48
## 5 PA14_00060  198.0543     -0.4497606 0.2497138  -1.8011040 7.168650e-02
## 6 PA14_00070  111.0326      0.1408506 0.3099986   0.4543588 6.495706e-01
##           padj     gene_id   seqnames start  end width strand    source type
## 1 6.090526e-47 gene1650835 chromosome   483 2027  1545      + PseudoCAP gene
## 2 4.006038e-18 gene1650837 chromosome  2056 3159  1104      + PseudoCAP gene
## 3 1.833081e-01 gene1650839 chromosome  3169 4278  1110      + PseudoCAP gene
## 4 4.831800e-47 gene1650841 chromosome  4275 6695  2421      + PseudoCAP gene
## 5 1.021604e-01 gene1650843 chromosome  7018 7791   774      - PseudoCAP gene
## 6 7.074361e-01 gene1650845 chromosome  7803 8339   537      - PseudoCAP gene
##   score phase Name         Dbxref       name Parent locus SK061_rpkm SK062_rpkm
## 1    NA     0      GeneID:4384099       dnaA         <NA>   7.915314   8.117840
## 2    NA     0      GeneID:4384100       dnaN         <NA>   8.403091   8.328621
## 3    NA     0      GeneID:4384101       recF         <NA>   6.919760   6.657980
## 4    NA     0      GeneID:4384102       gyrB         <NA>   7.947574   7.955889
## 5    NA     0      GeneID:4385186 PA14_00060         <NA>   6.474956   6.361612
## 6    NA     0      GeneID:4385187 PA14_00070         <NA>   5.867393   5.679355
##   SK063_rpkm SK082_rpkm SK083_rpkm SK084_rpkm SK061_cpm SK062_cpm SK063_cpm
## 1   8.104660   5.091606   4.892971   5.057935  8.540811  8.743613  8.730417
## 2   8.618680   6.719347   6.410163   7.087401  8.545430  8.470939  8.761074
## 3   7.126030   6.773662   6.474446   6.636371  7.069138  6.807123  7.275566
## 4   7.935406   5.265142   5.150084   4.930738  9.219743  9.228078  9.207545
## 5   6.462976   5.917922   5.957298   6.245160  6.110090  5.997131  6.098149
## 6   5.967636   6.328308   5.620217   6.119653  4.991538  4.806420  5.090371
##   SK082_cpm SK083_cpm SK084_cpm SK061_rawcounts SK062_rawcounts SK063_rawcounts
## 1  5.704210  5.503347  5.670183            1555            1680            1858
## 2  6.860797  6.551304  7.229142            1560            1390            1898
## 3  6.922914  6.623397  6.785493             558             436             675
## 4  6.518556  6.401640  6.178307            2492            2352            2588
## 5  5.555278  5.594467  5.881108             285             247             296
## 6  5.446699  4.748281  5.240425             129             106             145
##   SK082_rawcounts SK083_rawcounts SK084_rawcounts locusId   accession        GI
## 1             150             120             137 2194572 YP_788156.1 116053721
## 2             338             251             409 2194573 YP_788157.1 116053722
## 3             353             264             300 2194574 YP_788158.1 116053723
## 4             266             226             196 2194575 YP_788159.1 116053724
## 5             135             128             159 2194576 YP_788160.1 116053725
## 6             125              70             101 2194577 YP_788161.1 116053726
##   scaffoldId stop                                                  desc     COG
## 1       4582 2027 chromosomal replication initiator protein DnaA (NCBI)  COG593
## 2       4582 3159                 DNA polymerase III, beta chain (NCBI)  COG592
## 3       4582 4278        DNA replication and repair protein RecF (NCBI) COG1195
## 4       4582 6695                           DNA gyrase subunit B (NCBI)  COG187
## 5       4582 7018                       putative acyltransferase (NCBI)  COG204
## 6       4582 7803                putative histidinol-phosphatase (NCBI)  COG241
##   COGFun
## 1      L
## 2      L
## 3      L
## 4      L
## 5      I
## 6      E
##                                                                    COGDesc
## 1                            ATPase involved in DNA replication initiation
## 2                      DNA polymerase sliding clamp subunit (PCNA homolog)
## 3                         Recombinational DNA repair ATPase (RecF pathway)
## 4 Type IIA topoisomerase (DNA gyrase/topo II, topoisomerase IV), B subunit
## 5                           1-acyl-sn-glycerol-3-phosphate acyltransferase
## 6                          Histidinol phosphatase and related phosphatases
##                                                                                  TIGRFam
## 1                        TIGR00362 chromosomal replication initiator protein DnaA [dnaA]
## 2                                      TIGR00663 DNA polymerase III, beta subunit [dnaN]
## 3                               TIGR00611 DNA replication and repair protein RecF [recF]
## 4                                                 TIGR01059 DNA gyrase, B subunit [gyrB]
## 5                                                                                       
## 6 TIGR01656 histidinol-phosphate phosphatase domain,TIGR01662 HAD hydrolase, family IIIA
##                                                   TIGRRoles
## 1 DNA metabolism:DNA replication, recombination, and repair
## 2 DNA metabolism:DNA replication, recombination, and repair
## 3 DNA metabolism:DNA replication, recombination, and repair
## 4 DNA metabolism:DNA replication, recombination, and repair
## 5                                                          
## 6           Unknown function:Enzymes of unknown specificity
##                                                                                                                                                                                           GO
## 1                                                                                                                                     GO:0006270,GO:0006275,GO:0003688,GO:0017111,GO:0005524
## 2 GO:0006260,GO:0003677,GO:0003893,GO:0008408,GO:0016449,GO:0019984,GO:0003889,GO:0003894,GO:0015999,GO:0016450,GO:0003890,GO:0003895,GO:0016000,GO:0016451,GO:0003891,GO:0016448,GO:0016452
## 3                                                                                                                                     GO:0006281,GO:0005694,GO:0005524,GO:0017111,GO:0003697
## 4                                                                                                                                     GO:0006304,GO:0006265,GO:0005694,GO:0003918,GO:0005524
## 5                                                                                                                                                                      GO:0008152,GO:0003841
## 6                                                                                                                                                                      GO:0000105,GO:0004401
##         EC                                        ECDesc
## 1                                                       
## 2  2.7.7.7                  DNA-directed DNA polymerase.
## 3                                                       
## 4 5.99.1.3          DNA topoisomerase (ATP-hydrolyzing).
## 5 2.3.1.51 1-acylglycerol-3-phosphate O-acyltransferase.
## 6  3.1.3.-
write.csv(resdf, paste0(Sys.Date() , "_results_",exp, "_", treatment,".vs.", control,".csv"), row.names = FALSE)

9 Save DEG list

# filter for significantly differentially expressed genes (LFC >= 3, padj < pval_cutoff)
dim(resdf)
## [1] 5739   54
head(resdf)
##        Alias  baseMean log2FoldChange     lfcSE        stat       pvalue
## 1 PA14_00010  786.5229     -3.1219412 0.2140327 -14.5862786 3.434268e-48
## 2 PA14_00020  865.7989     -1.7706474 0.1993684  -8.8812834 6.609336e-19
## 3 PA14_00030  411.3546     -0.3275486 0.2205452  -1.4851766 1.374970e-01
## 4 PA14_00050 1164.6519     -2.9177741 0.1998135 -14.6024908 2.707692e-48
## 5 PA14_00060  198.0543     -0.4497606 0.2497138  -1.8011040 7.168650e-02
## 6 PA14_00070  111.0326      0.1408506 0.3099986   0.4543588 6.495706e-01
##           padj     gene_id   seqnames start  end width strand    source type
## 1 6.090526e-47 gene1650835 chromosome   483 2027  1545      + PseudoCAP gene
## 2 4.006038e-18 gene1650837 chromosome  2056 3159  1104      + PseudoCAP gene
## 3 1.833081e-01 gene1650839 chromosome  3169 4278  1110      + PseudoCAP gene
## 4 4.831800e-47 gene1650841 chromosome  4275 6695  2421      + PseudoCAP gene
## 5 1.021604e-01 gene1650843 chromosome  7018 7791   774      - PseudoCAP gene
## 6 7.074361e-01 gene1650845 chromosome  7803 8339   537      - PseudoCAP gene
##   score phase Name         Dbxref       name Parent locus SK061_rpkm SK062_rpkm
## 1    NA     0      GeneID:4384099       dnaA         <NA>   7.915314   8.117840
## 2    NA     0      GeneID:4384100       dnaN         <NA>   8.403091   8.328621
## 3    NA     0      GeneID:4384101       recF         <NA>   6.919760   6.657980
## 4    NA     0      GeneID:4384102       gyrB         <NA>   7.947574   7.955889
## 5    NA     0      GeneID:4385186 PA14_00060         <NA>   6.474956   6.361612
## 6    NA     0      GeneID:4385187 PA14_00070         <NA>   5.867393   5.679355
##   SK063_rpkm SK082_rpkm SK083_rpkm SK084_rpkm SK061_cpm SK062_cpm SK063_cpm
## 1   8.104660   5.091606   4.892971   5.057935  8.540811  8.743613  8.730417
## 2   8.618680   6.719347   6.410163   7.087401  8.545430  8.470939  8.761074
## 3   7.126030   6.773662   6.474446   6.636371  7.069138  6.807123  7.275566
## 4   7.935406   5.265142   5.150084   4.930738  9.219743  9.228078  9.207545
## 5   6.462976   5.917922   5.957298   6.245160  6.110090  5.997131  6.098149
## 6   5.967636   6.328308   5.620217   6.119653  4.991538  4.806420  5.090371
##   SK082_cpm SK083_cpm SK084_cpm SK061_rawcounts SK062_rawcounts SK063_rawcounts
## 1  5.704210  5.503347  5.670183            1555            1680            1858
## 2  6.860797  6.551304  7.229142            1560            1390            1898
## 3  6.922914  6.623397  6.785493             558             436             675
## 4  6.518556  6.401640  6.178307            2492            2352            2588
## 5  5.555278  5.594467  5.881108             285             247             296
## 6  5.446699  4.748281  5.240425             129             106             145
##   SK082_rawcounts SK083_rawcounts SK084_rawcounts locusId   accession        GI
## 1             150             120             137 2194572 YP_788156.1 116053721
## 2             338             251             409 2194573 YP_788157.1 116053722
## 3             353             264             300 2194574 YP_788158.1 116053723
## 4             266             226             196 2194575 YP_788159.1 116053724
## 5             135             128             159 2194576 YP_788160.1 116053725
## 6             125              70             101 2194577 YP_788161.1 116053726
##   scaffoldId stop                                                  desc     COG
## 1       4582 2027 chromosomal replication initiator protein DnaA (NCBI)  COG593
## 2       4582 3159                 DNA polymerase III, beta chain (NCBI)  COG592
## 3       4582 4278        DNA replication and repair protein RecF (NCBI) COG1195
## 4       4582 6695                           DNA gyrase subunit B (NCBI)  COG187
## 5       4582 7018                       putative acyltransferase (NCBI)  COG204
## 6       4582 7803                putative histidinol-phosphatase (NCBI)  COG241
##   COGFun
## 1      L
## 2      L
## 3      L
## 4      L
## 5      I
## 6      E
##                                                                    COGDesc
## 1                            ATPase involved in DNA replication initiation
## 2                      DNA polymerase sliding clamp subunit (PCNA homolog)
## 3                         Recombinational DNA repair ATPase (RecF pathway)
## 4 Type IIA topoisomerase (DNA gyrase/topo II, topoisomerase IV), B subunit
## 5                           1-acyl-sn-glycerol-3-phosphate acyltransferase
## 6                          Histidinol phosphatase and related phosphatases
##                                                                                  TIGRFam
## 1                        TIGR00362 chromosomal replication initiator protein DnaA [dnaA]
## 2                                      TIGR00663 DNA polymerase III, beta subunit [dnaN]
## 3                               TIGR00611 DNA replication and repair protein RecF [recF]
## 4                                                 TIGR01059 DNA gyrase, B subunit [gyrB]
## 5                                                                                       
## 6 TIGR01656 histidinol-phosphate phosphatase domain,TIGR01662 HAD hydrolase, family IIIA
##                                                   TIGRRoles
## 1 DNA metabolism:DNA replication, recombination, and repair
## 2 DNA metabolism:DNA replication, recombination, and repair
## 3 DNA metabolism:DNA replication, recombination, and repair
## 4 DNA metabolism:DNA replication, recombination, and repair
## 5                                                          
## 6           Unknown function:Enzymes of unknown specificity
##                                                                                                                                                                                           GO
## 1                                                                                                                                     GO:0006270,GO:0006275,GO:0003688,GO:0017111,GO:0005524
## 2 GO:0006260,GO:0003677,GO:0003893,GO:0008408,GO:0016449,GO:0019984,GO:0003889,GO:0003894,GO:0015999,GO:0016450,GO:0003890,GO:0003895,GO:0016000,GO:0016451,GO:0003891,GO:0016448,GO:0016452
## 3                                                                                                                                     GO:0006281,GO:0005694,GO:0005524,GO:0017111,GO:0003697
## 4                                                                                                                                     GO:0006304,GO:0006265,GO:0005694,GO:0003918,GO:0005524
## 5                                                                                                                                                                      GO:0008152,GO:0003841
## 6                                                                                                                                                                      GO:0000105,GO:0004401
##         EC                                        ECDesc
## 1                                                       
## 2  2.7.7.7                  DNA-directed DNA polymerase.
## 3                                                       
## 4 5.99.1.3          DNA topoisomerase (ATP-hydrolyzing).
## 5 2.3.1.51 1-acylglycerol-3-phosphate O-acyltransferase.
## 6  3.1.3.-
summary(is.na(resdf$log2FoldChange))
##    Mode   FALSE 
## logical    5739
summary(is.na(resdf$padj))
##    Mode   FALSE 
## logical    5739
sig_genes <- resdf[abs(resdf$log2FoldChange) >= logfc_cutoff & resdf$padj < pval_cutoff, ]  %>% drop_na(log2FoldChange)
sig_genes$direction <- ifelse(sig_genes$log2FoldChange > 0, "up", "down")

summary(is.na(sig_genes$log2FoldChange))
##    Mode   FALSE 
## logical    1403
summary(is.na(sig_genes$padj))
##    Mode   FALSE 
## logical    1403
dim(sig_genes)
## [1] 1403   55
head(sig_genes)
##        Alias   baseMean log2FoldChange     lfcSE       stat       pvalue
## 1 PA14_00010  786.52289      -3.121941 0.2140327 -14.586279 3.434268e-48
## 2 PA14_00050 1164.65188      -2.917774 0.1998135 -14.602491 2.707692e-48
## 3 PA14_00100   58.35071      -2.972299 0.4484006  -6.628669 3.387267e-11
## 4 PA14_00150   31.32745      -2.084440 0.5746102  -3.627573 2.860984e-04
## 5 PA14_00520   34.61451      -2.471521 0.5871134  -4.209615 2.558065e-05
## 6 PA14_00570   32.69422      -2.843659 0.5617334  -5.062293 4.142448e-07
##           padj     gene_id   seqnames start   end width strand    source type
## 1 6.090526e-47 gene1650835 chromosome   483  2027  1545      + PseudoCAP gene
## 2 4.831800e-47 gene1650841 chromosome  4275  6695  2421      + PseudoCAP gene
## 3 1.343219e-10 gene1650851 chromosome 12488 13435   948      - PseudoCAP gene
## 4 6.138616e-04 gene1650861 chromosome 16336 16608   273      - PseudoCAP gene
## 5 6.215070e-05 gene1650919 chromosome 54577 54699   123      + PseudoCAP gene
## 6 1.188343e-06 gene1650927 chromosome 59018 59704   687      + PseudoCAP gene
##   score phase Name         Dbxref       name Parent locus SK061_rpkm SK062_rpkm
## 1    NA     0      GeneID:4384099       dnaA         <NA>   7.915314   8.117840
## 2    NA     0      GeneID:4384102       gyrB         <NA>   7.947574   7.955889
## 3    NA     0      GeneID:4385190       glyQ         <NA>   5.210641   4.698845
## 4    NA     0      GeneID:4385195 PA14_00150         <NA>   5.643763   5.544972
## 5    NA     0      GeneID:4384761 PA14_00520         <NA>   7.157468   7.037571
## 6    NA     0      GeneID:4381739 PA14_00570         <NA>   4.758961   4.771188
##   SK063_rpkm SK082_rpkm SK083_rpkm SK084_rpkm SK061_cpm SK062_cpm SK063_cpm
## 1   8.104660   5.091606   4.892971   5.057935  8.540811  8.743613  8.730417
## 2   7.935406   5.265142   5.150084   4.930738  9.219743  9.228078  9.207545
## 3   5.135895   2.410146   2.505718   2.156835  5.135735  4.624848  5.061103
## 4   6.047883   4.390936   3.862043   2.939927  3.845603  3.751973  4.231788
## 5   7.031107   5.887614   4.461916   3.305598  4.204511  4.090548  4.084417
## 6   4.409378   2.161430   2.080155   2.236586  4.241417  4.253442  3.898365
##   SK082_cpm SK083_cpm SK084_cpm SK061_rawcounts SK062_rawcounts SK063_rawcounts
## 1  5.704210  5.503347  5.670183            1555            1680            1858
## 2  6.518556  6.401640  6.178307            2492            2352            2588
## 3  2.347917  2.442544  2.097431             143              93             142
## 4  2.690309  2.231636  1.496680              56              49              78
## 5  3.028401  1.843045  1.065680              73              63              70
## 6  1.759731  1.686163  1.828102              75              71              61
##   SK082_rawcounts SK083_rawcounts SK084_rawcounts locusId   accession        GI
## 1             150             120             137 2194572 YP_788156.1 116053721
## 2             266             226             196 2194575 YP_788159.1 116053724
## 3              12              12               9 2194580 YP_788164.1 116053729
## 4              16              10               5 2194585 YP_788169.1 116053734
## 5              21               7               3 2194614 YP_788198.1 116053763
## 6               7               6               7 2194618 YP_788202.1 116053767
##   scaffoldId  stop                                                  desc
## 1       4582  2027 chromosomal replication initiator protein DnaA (NCBI)
## 2       4582  6695                           DNA gyrase subunit B (NCBI)
## 3       4582 12488             glycyl-tRNA synthetase alpha chain (NCBI)
## 4       4582 16336                           hypothetical protein (NCBI)
## 5       4582 54699                           hypothetical protein (NCBI)
## 6       4582 59704                           putative lipoprotein (NCBI)
##       COG COGFun
## 1  COG593      L
## 2  COG187      L
## 3  COG752      J
## 4               
## 5               
## 6 COG1462      M
##                                                                    COGDesc
## 1                            ATPase involved in DNA replication initiation
## 2 Type IIA topoisomerase (DNA gyrase/topo II, topoisomerase IV), B subunit
## 3                                    Glycyl-tRNA synthetase, alpha subunit
## 4                                                                         
## 5                                                                         
## 6          Uncharacterized protein involved in formation of curli polymers
##                                                           TIGRFam
## 1 TIGR00362 chromosomal replication initiator protein DnaA [dnaA]
## 2                          TIGR01059 DNA gyrase, B subunit [gyrB]
## 3            TIGR00388 glycine--tRNA ligase, alpha subunit [glyQ]
## 4                                                                
## 5                                                                
## 6                                                                
##                                                   TIGRRoles
## 1 DNA metabolism:DNA replication, recombination, and repair
## 2 DNA metabolism:DNA replication, recombination, and repair
## 3                     Protein synthesis:tRNA aminoacylation
## 4                                                          
## 5                                                          
## 6                                                          
##                                                       GO       EC
## 1 GO:0006270,GO:0006275,GO:0003688,GO:0017111,GO:0005524         
## 2 GO:0006304,GO:0006265,GO:0005694,GO:0003918,GO:0005524 5.99.1.3
## 3            GO:0006426,GO:0009345,GO:0005524,GO:0004820 6.1.1.14
## 4                                                                
## 5                                                                
## 6                                  GO:0015031,GO:0042597         
##                                 ECDesc direction
## 1                                           down
## 2 DNA topoisomerase (ATP-hydrolyzing).      down
## 3                Glycine--tRNA ligase.      down
## 4                                           down
## 5                                           down
## 6                                           down
write.csv(sig_genes, paste0(Sys.Date(),"_significantgenes_", exp, "_", treatment,".vs.", control,".csv"), row.names = FALSE)

10 PCA plot

vst <- DESeqDataSetFromMatrix(countData = rawcounts, colData = metadata, design = ~ exp_11 + batch)
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
# Variance stabilizing transformations
vsd <- vst(vst, blind = TRUE)

# PCA plot
pcaData <- plotPCA(vsd, intgroup = c("exp_11", "batch", "condition"), returnData = TRUE)
## using ntop=500 top features by variance
percentVar <- round(100 * attr(pcaData, "percentVar"))

# Calculate the maximum absolute value among PC1 and PC2
max_abs_val <- max(c(abs(min(pcaData$PC1)), abs(min(pcaData$PC2)), abs(max(pcaData$PC1)), abs(max(pcaData$PC2))))

axis_size <- 10
theme_set(theme_minimal(base_family = "Arial"))

pca <- ggplot(pcaData, aes(PC1, PC2, color = condition, shape = condition)) +
  geom_point(size = 3) +
  xlab(paste0("PC1: ", percentVar[1], "% variance")) +
  ylab(paste0("PC2: ", percentVar[2], "% variance")) +
  scale_colour_brewer(palette = "Dark2") +
  scale_shape_manual(values = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16)) +
  coord_fixed(ratio = 1, xlim = c(-max_abs_val, max_abs_val), ylim = c(-max_abs_val, max_abs_val)) +
  labs(color = "condition", shape = "condition") +
  guides(shape = guide_legend(title = "Condition"),  # Adjust size of shape legend
    color = guide_legend(title = "Condition")) +
  theme(
    panel.grid = element_blank(),
    panel.background = element_blank(),
    axis.text = element_text(size = axis_size),
    axis.ticks = element_line(linewidth = 0.2),
     axis.line = element_line(linewidth = 0, color = "black"),  # Set the axis line linewidth
   panel.border = element_rect(fill = NA, linewidth = 0.5),  # Set the panel border linewidth
  legend.key = element_rect(fill = "white", color = NA)  # Set legend key background to white and remove border
  )

print(pca)

cairo_pdf(paste0(Sys.Date(),"_pca_",exp, "_", treatment,".vs.", control,".pdf", sep = ""))
pca
dev.off()
## quartz_off_screen 
##                 2

11 Volcano plot

sequence <- c("treatment")

for (i in sequence) {
    set.seed(0)
    res <- results(dds_final, contrast = c("exp_11", i, "control"))

  gp2 <- wes_palettes$Zissou1
  key <- read.csv("../../raw-data/pa14_gene_key.csv",  header = TRUE)

  rownames(key) <- key$gene_id

  res <- merge(key, as.data.frame(res), by=0)
  res <- as.data.frame(res)
  res$Alias <- rownames(res)

 # Define gene list
  alg_genes <- c(
    "algA", "algF", "algJ", "algI", "algL", "algX", "algG", "algE",
    "algK", "alg44", "alg8", "algD", "mucC", "mucB", "algU", "algW",
    "algP", "algQ", "algR", "algZ", "algC", "algB"
  )

logfc_cutoff
keyvals <- ifelse(
 res$log2FoldChange <= -logfc_cutoff & res$padj <= pval_cutoff, gp2[1],
  ifelse(res$log2FoldChange >= logfc_cutoff & res$padj <= pval_cutoff, gp2[5],
         'lightgrey'))
keyvals[is.na(keyvals)] <- 'lightgrey'
names(keyvals)[keyvals == gp2[5]] <- 'Up-regulated'
names(keyvals)[keyvals == 'lightgrey'] <- 'NS'
names(keyvals)[keyvals == gp2[1]] <- 'Down-regulated'

volcano <- EnhancedVolcano(res,
                     lab = res[["name"]],
                     selectLab = alg_genes,
                     x = 'log2FoldChange',
                     y = 'padj',
                     xlab = bquote(~Log[2]~ 'fold change'),
                     pCutoff = pval_cutoff,
                     FCcutoff = logfc_cutoff,
                     cutoffLineType = "dashed",
                     pointSize = 1,
                     labSize = 3.0,
                     axisLabSize = 14,
                     labCol = 'black',
                     labFace = 'bold',
                     #col = c('grey', gp2[2], gp2[3], gp2[1]),
                     colCustom = keyvals,
                     colAlpha = 5/5,
                     legendPosition = 'bottom',
                     legendLabSize = 12,
                     legendIconSize = 4.0,
                     drawConnectors = TRUE,
                     arrowheads = FALSE,
                     widthConnectors = 0.2,
                     gridlines.major = FALSE,
                     borderWidth = 0.3,
                     gridlines.minor = FALSE,
                     title = paste0(treatment," vs ", control, sep = ""),
                     subtitle = "",
                     #caption = paste0("total = ", nrow(toptable), " variables"),
                     caption = paste0("Volcano plot with log2 FC of ", logfc_cutoff," and P value of ",pval_cutoff, "."),
                     max.overlaps = 30,
                     #ylim = c(0,10),
                     colConnectors = 'black') + 
      theme(axis.ticks=element_line(linewidth =0.2)) # +
 # scale_x_continuous(breaks = c(-5, -2.5, 0, 2.5, 5), limits = c(-6,6))

print(volcano)

cairo_pdf(paste0(Sys.Date(),"_volcano_", exp, "_", treatment,".vs.", control,".pdf", sep = ""))
print(volcano)
dev.off()
}
## Warning: One or more p-values is 0. Converting to 10^-1 * current lowest
## non-zero p-value...

12 MA plot

sequence <- c("treatment")

for (i in sequence) {
    set.seed(0)

res <- results(dds_final, contrast = c("exp_11", i, "control"))

gp2 <- wes_palettes$Zissou1

res$gene_id <- rownames(res)
res <- as.data.frame(res)
head(res)
res <- full_join(key, as.data.frame(res), by="gene_id")
rownames(res) <- res$Alias

any(is.na(res))  # Check if there are any NAs

res <- res %>% filter(!is.na(baseMean))


ma <- ggmaplot(res, 
               fdr = pval_cutoff, 
               fc = 2^(logfc_cutoff), 
               size = 1,
               palette = c((gp2[5]), (gp2[1])),
               select.top.method = "fc",
             #  ylim = c(-15,15),
               legend = "bottom", 
               title = paste0(treatment," vs ", control, sep = ""),
              subtitle = paste0("Volcano plot with log2 FC of ", logfc_cutoff," and P value of ",pval_cutoff, "."),
               font.label = c( 12), 
             label.rectangle = FALSE,
              genenames = res$name,
              font.legend = c(12),
               font.main = c("bold", 14),
               alpha = 1.5,
               ggtheme = ggplot2::theme(panel.grid =element_blank(),
                                        axis.text = element_text(size = 12),
                                        axis.title = element_text(size = 14),
                                        panel.background = element_blank(),
                                        panel.border = element_rect(fill = NA))) 
print(ma)

cairo_pdf(paste0(Sys.Date(),"_ma_", exp, "_", treatment,".vs.", control,".pdf"))
print(ma)
dev.off()
}

LS0tCnRpdGxlOiAiMjAyNV9leHBfMTFfYWxnaW5hdGVfREVHIgphdXRob3I6ICJOb3VyIEVsIEh1c3NlaW5pIgpkYXRlOiAiMjAyNSIKb3V0cHV0OiAKICBodG1sX2RvY3VtZW50OgogICAgdG9jOiB5ZXMKICAgIHRvY19mbG9hdDogeWVzCiAgICB0b2NfY29sbGFwc2VkOiBmYWxzZSAgIyBFbnN1cmVzIFRPQyBpcyBleHBhbmRlZCBieSBkZWZhdWx0CiAgICBudW1iZXJfc2VjdGlvbnM6IHllcwogICAgY29kZV9mb2xkaW5nOiBzaG93CiAgICBjb2RlX2Rvd25sb2FkOiB0cnVlCiAgICBoaWdobGlnaHQ6IGthdGUKICAgIHRoZW1lOgogICAgICBib290c3dhdGNoOiB6ZXBoeXIKICAgIGNzczogc3R5bGUuY3NzICAjIExpbmsgdG8geW91ciBjdXN0b20gQ1NTIGZpbGUKZWRpdG9yX29wdGlvbnM6IAogIGNodW5rX291dHB1dF90eXBlOiBjb25zb2xlCiAgbWFya2Rvd246IAogICAgd3JhcDogNzIKLS0tCgpgYGB7ciBsb2FkX2xpYnJhcmllcywgZWNobz1GQUxTRSwgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KbGlicmFyeSgiREVTZXEyIikKbGlicmFyeSgiZ2dwbG90MiIpCmxpYnJhcnkoImthYmxlRXh0cmEiKQpsaWJyYXJ5KCJkcGx5ciIpCmxpYnJhcnkoImdnZm9yY2UiKQpsaWJyYXJ5KCJSQ29sb3JCcmV3ZXIiKQpsaWJyYXJ5KCJwaGVhdG1hcCIpCmxpYnJhcnkoIndlc2FuZGVyc29uIikKbGlicmFyeSgiRW5oYW5jZWRWb2xjYW5vIikKbGlicmFyeSgiZ2dwdWJyIikKbGlicmFyeSgiaHBnbHRvb2xzIikKbGlicmFyeSgiZHBseXIiKQpsaWJyYXJ5KCJ0aWJibGUiKSAKbGlicmFyeSgiZWRnZVIiKQpsaWJyYXJ5KCJnZW5lZmlsdGVyIikKbGlicmFyeSgicGhlYXRtYXAiKQpsaWJyYXJ5KCJkcGx5ciIpCmxpYnJhcnkoInRpZHl2ZXJzZSIpCmxpYnJhcnkoIkRFU2VxMiIpCmxpYnJhcnkoIkNvbXBsZXhIZWF0bWFwIikKbGlicmFyeSgiZ29zZXEiKQpsaWJyYXJ5KCJwbHlyIikKbGlicmFyeSgiZ2VvbXRleHRwYXRoIikKYGBgCgojIE92ZXJ2aWV3CgpgYGB7ciBzYW1wbGVfdGFibGUsIGVjaG89RkFMU0V9CiMgTG9hZCBtZXRhZGF0YSBDU1YgZmlsZSBpbnRvIGEgZGF0YWZyYW1lLCB1c2luZyB0aGUgZmlyc3QgY29sdW1uIGFzIHJvdyBuYW1lcwptZXRhZGF0YSA8LSByZWFkLmNzdigiLi4vLi4vcmF3LWRhdGEvbWV0YWRhdGFfYWxnaW5hdGVfdnRsLmNzdiIsIHJvdy5uYW1lcyA9IDEsIGhlYWRlciA9IFRSVUUpCgojIEZpbHRlciB0aGUgZGF0YSB0byBpbmNsdWRlIG9ubHkgcm93cyB3aGVyZSBleHBfMTEgaXMgZWl0aGVyICJjb250cm9sIiBvciAidHJlYXRtZW50IgojIFRoZW4gcmVtb3ZlIG90aGVyIGV4cGVyaW1lbnRzIChleHBfMiB0byBleHBfNikKbWV0YWRhdGEgPC0gZmlsdGVyKG1ldGFkYXRhLCBleHBfMTEgJWluJSBjKCJjb250cm9sIiwgInRyZWF0bWVudCIpKSAlPiUKICAgIHNlbGVjdChiYXRjaCwgY29uZGl0aW9uLCBtZWRpYSwgdHJlYXRtZW50LCBpbmR1Y2VkLCBleHBfMTEpIAoKIyBEaXNwbGF5IHRoZSBjbGVhbmVkIG1ldGFkYXRhIHRhYmxlIHVzaW5nIGtuaXRyOjprYWJsZSBhbmQga2FibGVFeHRyYSBzdHlsaW5nCm1ldGFkYXRhICU+JQogIGtibCgpICU+JSAgIyBDcmVhdGUgYSBiYXNpYyB0YWJsZQogIGthYmxlX2NsYXNzaWMoZnVsbF93aWR0aCA9IFRSVUUpICAjIEFwcGx5IGEgY2xhc3NpYyB0YWJsZSBzdHlsZSB3aXRoIGZ1bGwgd2lkdGgKYGBgIAoKIyBQcm9qZWN0IHZhcmlhYmxlcwpgYGB7ciBwcm9qZWN0X3BhcmFtZXRlcnMsIGVjaG89RkFMU0V9CiMgRGVmaW5lIHRoZSBsYWJlbHMgZm9yIGNvbnRyb2wgYW5kIHRyZWF0bWVudCBncm91cHMgdXNlZCBpbiB0aGUgZXhwZXJpbWVudApoZWFkKG1ldGFkYXRhKQojIEV4dHJhY3QgY29udHJvbCBhbmQgdHJlYXRtZW50IGNvbmRpdGlvbnMgYmFzZWQgb24gYGV4cF8xMWAgbGFiZWxzCmNvbnRyb2wgPC0gbWV0YWRhdGFbbWV0YWRhdGEkZXhwXzExID09ICJjb250cm9sIiwgXQp0cmVhdG1lbnQgPC0gbWV0YWRhdGFbbWV0YWRhdGEkZXhwXzExID09ICJ0cmVhdG1lbnQiLCBdCgpjb250cm9sIDwtIHVuaXF1ZShjb250cm9sJGNvbmRpdGlvbikKdHJlYXRtZW50IDwtIHVuaXF1ZSh0cmVhdG1lbnQkY29uZGl0aW9uKQoKcHJpbnQocGFzdGUwKCJUaGUgY29udHJvbCBpczogIiwgY29udHJvbCkpCnByaW50KHBhc3RlMCgiVGhlIHRyZWF0bWVudCBpczogIiwgdHJlYXRtZW50KSkKCiMgU2V0IGRlZmF1bHQgdGV4dCBzaXplcyBmb3IgcGxvdHMKdGV4dC5zaXplIDwtIDEyICAgIApheGlzLnRleHQuc2l6ZSA8LSAxNAoKbG9nZmNfY3V0b2ZmIDwtIDIKcHZhbF9jdXRvZmYgPC0gMC4wMQoKIyBTZXQgdmFyaWFibGVzIGZvciB0aGUgcHJvamVjdCB1c2VkIGluIGZpbGVuYW1lcyBhbmQgcGxvdCB0aXRsZXMKcHJvamVjdCA8LSAib3JuX2FsZ2luYXRlX3JuYXNlcSIKZXhwIDwtICJleHBfMTEiICAgICAgICAgICAgICAgICAgICAKYGBgCgpgYGB7ciByYXdkYXRhfQojIExvYWQgcmF3IGNvdW50IGRhdGEKcmF3Y291bnRzIDwtIHJlYWQuY3N2KCIuLi8uLi9yYXctZGF0YS9odHNlcV9TS191bmZpbHRlcmVkX3JldmVyc2Vfc3RyYW5kX2NvdW50LmNzdiIsIAogICAgICAgICAgICAgICAgICAgICAgcm93Lm5hbWVzID0gMSwgaGVhZGVyID0gVFJVRSkKCiMgQXBwZW5kIGEgcHJlZml4IHRvIHRoZSBiYXRjaCBhbmQgcnVuIGNvbHVtbnMgdG8gdHVybiB0byBjaGFyYWN0ZXIgY2xhc3MKbWV0YWRhdGFbWyJiYXRjaCJdXSA8LSBwYXN0ZTAoImIiLCBtZXRhZGF0YVtbImJhdGNoIl1dKQptZXRhZGF0YVtbInJ1biJdXSA8LSBwYXN0ZTAoInIiLCBtZXRhZGF0YVtbInJ1biJdXSkKCiMgU3Vic2V0IHRoZSByYXdjb3VudHMgZGF0YWZyYW1lIHRvIG9ubHkgaW5jbHVkZSBzYW1wbGVzIHRoYXQgYXJlIHByZXNlbnQgaW4gZXhwZXJpbWVudApzYW1wbGVzIDwtIHJvd25hbWVzKG1ldGFkYXRhKQpyYXdjb3VudHMgPC0gc2VsZWN0KHJhd2NvdW50cywgYWxsX29mKHNhbXBsZXMpKQoKIyBFbnN1cmUgdGhhdCB0aGUgb3JkZXIgb2Ygc2FtcGxlcyBpbiBtZXRhZGF0YSBtYXRjaGVzIHRoZSBvcmRlciBvZiBjb2x1bW5zIGluIHJhd2NvdW50cwphbGwocm93bmFtZXMobWV0YWRhdGEpID09IGNvbG5hbWVzKHJhd2NvdW50cykpICAjIFNob3VsZCByZXR1cm4gVFJVRSBpZiBvcmRlciBpcyBjb3JyZWN0CgojIERlZmluZSB0aGUgZGVzaWduIGZvcm11bGEgZm9yIERFU2VxMiwgaW5jbHVkaW5nIG1haW4gY29uZGl0aW9uIGFuZCBiYXRjaCBlZmZlY3QKZGVzaWduX2Zvcm11bGEgPC0gfiBleHBfMTEgKyBiYXRjaAoKIyBDcmVhdGUgREVTZXEyIGRhdGFzZXQgb2JqZWN0IGZyb20gdGhlIGZpbHRlcmVkIGNvdW50IG1hdHJpeCBhbmQgbWV0YWRhdGEKZGRzIDwtIERFU2VxRGF0YVNldEZyb21NYXRyaXgoCiAgY291bnREYXRhID0gcmF3Y291bnRzLAogIGNvbERhdGEgPSBtZXRhZGF0YSwKICBkZXNpZ24gPSBkZXNpZ25fZm9ybXVsYQopCmBgYAoKCiMgTGlicmFyeSBzaXplIHBlciBzYW1wbGUKYGBge3IgbGlicmFyeV9zaXplfQojIFN1bSB1cCB0aGUgY291bnRzIGZvciBlYWNoIHNhbXBsZQpyYXdjb3VudHNfc3VtcyA8LSBjb2xTdW1zKHJhd2NvdW50cykKCiMgVHVybiBpdCBpbnRvIGEgZGF0YSBmcmFtZQpkYXRhIDwtIGRhdGEuZnJhbWUoCiAgQ2F0ZWdvcnkgPSBuYW1lcyhyYXdjb3VudHNfc3VtcyksCiAgQ291bnQgPSBhcy5udW1lcmljKHJhd2NvdW50c19zdW1zKQopCgojIE1ha2Ugc3VyZSBzYW1wbGUgb3JkZXIgaXMgY29uc2lzdGVudApkYXRhIDwtIGRhdGEgJT4lCiAgbXV0YXRlKENhdGVnb3J5ID0gZmFjdG9yKENhdGVnb3J5LCBsZXZlbHMgPSBzb3J0KHVuaXF1ZShDYXRlZ29yeSkpKSkKCiMgUGxvdCBiYXIgY2hhcnQgb2YgbGlicmFyeSBzaXplcwpwIDwtIGdncGxvdChkYXRhLCBhZXMoeCA9IENhdGVnb3J5LCB5ID0gQ291bnQsIGZpbGwgPSBDYXRlZ29yeSkpICsKICBnZW9tX2JhcihzdGF0ID0gImlkZW50aXR5IikgKwogICMgY29vcmRfZmxpcCgpICMgZmxpcCBpZiB0b28gbWFueSBzYW1wbGVzIHRvIHJlYWQgbGFiZWxzIGVhc2lseQogIGxhYnMoeSA9ICJMaWJyYXJ5IHNpemUgKE0pIiwgeCA9IE5VTEwpICsgICMgcmVtb3ZlIHggbGFiZWwsIHJlbmFtZSB5IGF4aXMKICAjIHNjYWxlX2ZpbGxfYnJld2VyKHBhbGV0dGUgPSAiRGFyazIiKSArICAjIG9wdGlvbmFsIGNvbG9yIHBhbGV0dGUKICBzY2FsZV95X2NvbnRpbnVvdXMoCiAgICBsYWJlbHMgPSBzY2FsZXM6OmxhYmVsX251bWJlcihzY2FsZSA9IDFlLTYpLCAgIyBzaG93IGNvdW50cyBpbiBtaWxsaW9ucwogICAgbGltaXRzID0gYygwLCAobWF4KGRhdGEkQ291bnQpICsgMWU2KSkgICAgICAgICMgYWRkIGEgYml0IG9mIHNwYWNlIG9uIHRvcAogICkgKwogIHRoZW1lKAogICAgcGFuZWwuZ3JpZCA9IGVsZW1lbnRfYmxhbmsoKSwKICAgIHBhbmVsLmJhY2tncm91bmQgPSBlbGVtZW50X2JsYW5rKCksCiAgICBheGlzLnRleHQgPSBlbGVtZW50X3RleHQoc2l6ZSA9IHRleHQuc2l6ZSwgYW5nbGUgPSA0NSwgaGp1c3QgPSAxKSwgICMgc2xhbnQgbGFiZWxzCiAgICBheGlzLnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSBheGlzLnRleHQuc2l6ZSksCiAgICBsZWdlbmQudGV4dCA9IGVsZW1lbnRfdGV4dChzaXplID0gdGV4dC5zaXplKSwKICAgIGxlZ2VuZC50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gYXhpcy50ZXh0LnNpemUpLAogICAgcGFuZWwuYm9yZGVyID0gZWxlbWVudF9yZWN0KGZpbGwgPSBOQSksCiAgICBsZWdlbmQucG9zaXRpb24gPSAibm9uZSIsICAjIG5vIGxlZ2VuZCBuZWVkZWQgc2luY2UgZWFjaCBiYXIgaXMgYSBkaWZmZXJlbnQgc2FtcGxlCiAgICBheGlzLnRpY2tzID0gZWxlbWVudF9saW5lKGxpbmV3aWR0aCA9IDAuMiksCiAgICBsZWdlbmQua2V5ID0gZWxlbWVudF9ibGFuaygpCiAgKQoKIyBTaG93IHRoZSBwbG90CnAKCiMgU2F2ZSB0aGUgcGxvdCBhcyBhIFBERiB3aXRoIHRvZGF54oCZcyBkYXRlIGFuZCBwcm9qZWN0IG5hbWUKZ2dzYXZlKHBhc3RlMChTeXMuRGF0ZSgpLCAiX2xpYnJhcnlfc2l6ZV8iLCBleHAsICJfIiwgdHJlYXRtZW50LCIudnMuIiwgY29udHJvbCwgIi5wZGYiKSwgZGV2aWNlID0gInBkZiIpCmBgYAoKIyBFeHByZXNzaW9uIGRpc3RyaWJ1dGlvbgoKYGBge3IgZXhwcmVzc2lvbl9kaXN0cmlidXRpb259CgojIFN0YWNrIHJhd2NvdW50cyBpbnRvIHNpbmdsZSB2ZWN0b3IgCmRhdCA8LSBzdGFjayhyYXdjb3VudHMpCgojIEZvcm1hdCBleHByZXNzaW9uIGRpc3RyaWJ1dGlvbiBwbG90IHdpdGggdGhlIGV4cHJlc3Npb24gbGV2ZWwgYW5kIGRlbnNpdHkgCmV4cHJlc3Npb25fZGlzdCA8LSBnZ3Bsb3QoZGF0LCBhZXMoeCA9IGxvZzIodmFsdWVzKzEpLCAgZmlsbD1pbmQsIGxhYmVsID0gaW5kICkpICsgCiAgZ2VvbV9kZW5zaXR5KGFscGhhID0gMC41LCBwb3NpdGlvbiA9ICJpZGVudGl0eSIpICsgCiAgeWxhYigiRGVuc2l0eSIpICsgCiAgZ2VvbV90ZXh0ZGVuc2l0eShoanVzdCA9ICJ5bWF4Iiwgdmp1c3QgPSAwLAogICAgICAgICAgIHRleHRfb25seSA9IFRSVUUsIHRleHRfc21vb3RoaW5nID0gMjApICsgCiAgICB0aGVtZShwYW5lbC5ncmlkID0gZWxlbWVudF9ibGFuaygpLAogICAgICAgIHBhbmVsLmJhY2tncm91bmQgPSBlbGVtZW50X2JsYW5rKCksCiAgICAgICAgYXhpcy50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSB0ZXh0LnNpemUpLAogICAgICAgIGF4aXMudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IGF4aXMudGV4dC5zaXplKSwKICAgICAgICBsZWdlbmQudGV4dCA9IGVsZW1lbnRfdGV4dChzaXplID0gdGV4dC5zaXplKSwKICAgICAgICBsZWdlbmQudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IGF4aXMudGV4dC5zaXplKSwKICAgICAgICBwYW5lbC5ib3JkZXIgPSBlbGVtZW50X3JlY3QoZmlsbCA9IE5BKSwKICAgICAgICBsZWdlbmQucG9zaXRpb24gPSAibm9uZSIsCiAgICAgICAgYXhpcy50aWNrcyA9IGVsZW1lbnRfbGluZShsaW5ld2lkdGggPSAwLjIpLAogICAgICAgIGxlZ2VuZC5rZXkgPSBlbGVtZW50X2JsYW5rKCkpICsKICBjb29yZF9jYXJ0ZXNpYW4oeGxpbSA9IGMoMCwgMTUpKQoKZXhwcmVzc2lvbl9kaXN0Cgojc2F2ZSBmaWd1cmUKY2Fpcm9fcGRmKHBhc3RlMChTeXMuRGF0ZSgpLCAiX2V4cHJlc3Npb25fZGlzdF8iLCBleHAsICJfIiwgdHJlYXRtZW50LCAiLnZzLiIsIGNvbnRyb2wsICIucGRmIikpCmV4cHJlc3Npb25fZGlzdApkZXYub2ZmKCkgCgpgYGAKCiMgTm9uLXplcm8gZ2VuZXMgb2JzZXJ2ZWQKCmBgYHtyIG5vbi16ZXJvX2dlbmVzfQojIFByZXBhcmUgZGF0YSBieSBnZXR0aW5nIGxpYnJhcnkgc2l6ZXMgYW5kIG5vbi16ZXJvIGNvdW50cwp4IDwtIGNvbFN1bXMocmF3Y291bnRzKQp5IDwtIGNvbFN1bXMocmF3Y291bnRzICE9IDApIApkZiA8LSBkYXRhLmZyYW1lKHgseSkKZGYkaW5kIDwtIHJvd25hbWVzKGRmKQoKIyBDYWxjdWxhdGUgYXhpcyBsaW1pdHMKeF9taW4gPC0gbWluKGRmJHgpCnhfbWF4IDwtIG1heChkZiR4KQp5X21pbiA8LSBtaW4oZGYkeSkKeV9tYXggPC0gbWF4KGRmJHkpCgojIE5vbi16ZXJvIHBsb3QKbm9uemVybyA8LSBnZ3Bsb3QoZGF0YSA9IGRmLCBhZXMoeCA9IHgsIHkgPSB5LCBjb2xvciA9IGluZCwgbGFiZWwgPSBpbmQpKSArCiAgZ2VvbV9wb2ludChzaXplID0gMykgKwogIGdlb21fdGV4dChoanVzdCA9IC0wLjEsIHZqdXN0ID0gLTAuMSwgY29sb3IgPSAiYmxhY2siLCBzaXplID0gNSkgKwogIHNjYWxlX3lfY29udGludW91cyhsaW1pdHM9Yyh5X21pbiwgeV9tYXgpKSArIAogIHNjYWxlX3hfY29udGludW91cygKICAgIGxpbWl0cyA9IGMoeF9taW4tMWUzLCB4X21heCs5ZTUpLAogICAgbGFiZWxzID0gc2NhbGVzOjpjb21tYV9mb3JtYXQoc2NhbGUgPSAxZS02LCAgYWNjdXJhY3kgPSAxKSAgIyBGb3JtYXQgbGFiZWxzIHdpdGhvdXQgZGVjaW1hbHMKICApICsKICB4bGFiKCJNaWxsaW9ucyBvZiByZWFkcyIpICsgeWxhYigiTnVtYmVyIG9mIG5vbi16ZXJvIGdlbmVzIG9ic2VydmVkIikgKwogIGxhYnMoY29sb3IgPSAiU2FtcGxlIikgKyAKICAgIHRoZW1lKHBhbmVsLmdyaWQgPSBlbGVtZW50X2JsYW5rKCksCiAgICAgICAgcGFuZWwuYmFja2dyb3VuZCA9IGVsZW1lbnRfYmxhbmsoKSwKICAgICAgICBheGlzLnRleHQgPSBlbGVtZW50X3RleHQoc2l6ZSA9IHRleHQuc2l6ZSksCiAgICAgICAgYXhpcy50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gYXhpcy50ZXh0LnNpemUpLAogICAgICAgIGxlZ2VuZC50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSB0ZXh0LnNpemUpLAogICAgICAgIGxlZ2VuZC50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gYXhpcy50ZXh0LnNpemUpLAogICAgICAgIHBhbmVsLmJvcmRlciA9IGVsZW1lbnRfcmVjdChmaWxsID0gTkEpLAogICAgICAgIGxlZ2VuZC5wb3NpdGlvbiA9ICJub25lIiwKICAgICAgICBheGlzLnRpY2tzID0gZWxlbWVudF9saW5lKGxpbmV3aWR0aCA9IDAuMiksCiAgICAgICAgbGVnZW5kLmtleSA9IGVsZW1lbnRfYmxhbmsoKSkgKwogIGNvb3JkX2NhcnRlc2lhbigpCgojIFNob3cgbm9uemVybyBwbG90Cm5vbnplcm8KCiMgU2F2ZSBwbG90CgpjYWlyb19wZGYocGFzdGUwKFN5cy5EYXRlKCksICJfbm9uLXplcm9fIiwgZXhwLCAiXyIsIHRyZWF0bWVudCwgIi52cy4iLCBjb250cm9sLCAiLnBkZiIpKQpub256ZXJvCmRldi5vZmYoKSAKCmBgYAoKIyBSUEtNIGFuZCBDUE0gdHJhbnNmb3JtYXRpb25zCmBgYHtyIGRhdGFfdHJhbnNmb3JtYXRpb25zfQojIExvYWQgR0ZGIGFubm90YXRpb24gCmdmZl9hbm5vdCA8LSBsb2FkX2dmZl9hbm5vdGF0aW9ucygiLi4vLi4vcmF3LWRhdGEvcGFlcnVnaW5vc2FfcGExNC5nZmYiLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGlkX2NvbCA9ICJnZW5lX2lkIiwgdHlwZSA9ICJnZW5lIikKcm93bmFtZXMoZ2ZmX2Fubm90KSA8LSBnZmZfYW5ub3QkZ2VuZV9pZAoKaGVhZChnZmZfYW5ub3QpCmhlYWQocmF3Y291bnRzKQpkaW0oZ2ZmX2Fubm90KQpkaW0ocmF3Y291bnRzKQoKc3VtbWFyeShyb3duYW1lcyhnZmZfYW5ub3QpID09IHJvd25hbWVzKHJhd2NvdW50cykpCmdmZl9hbm5vdCA8LSBnZmZfYW5ub3Rbcm93bmFtZXMocmF3Y291bnRzKSxdCnN1bW1hcnkocm93bmFtZXMoZ2ZmX2Fubm90KSA9PSByb3duYW1lcyhyYXdjb3VudHMpKQoKCiMgQ3JlYXRlIGFuIGV4cHJlc3Npb24gb2JqZWN0IHdpdGggY291bnRzLCBtZXRhZGF0YSwgYW5kIGdlbmUgaW5mbwpleHB0IDwtIGNyZWF0ZV9leHB0KGNvdW50X2RhdGFmcmFtZSA9IHJhd2NvdW50cywKICAgICAgICAgICAgICAgICAgICBtZXRhZGF0YSA9IG1ldGFkYXRhLCAKICAgICAgICAgICAgICAgICAgICBnZW5lX2luZm8gPSBnZmZfYW5ub3QpCgojIE5vcm1hbGl6ZSBjb3VudHMgbG9nMiBDUE0KY3BtX25vcm1hbGl6ZWRfZXhwIDwtIG5vcm1hbGl6ZV9leHB0KGV4cHQsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdHJhbnNmb3JtID0gImxvZzIiLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG5vcm0gPSAicmF3IiwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjb252ZXJ0ID0gImNwbSIsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgcm93X21pbiA9ICIxMCIpCgojIyMjIyMjIENhbGN1bGF0ZSBSUEtNICMjIyMjIyMjCgojIEdyYWIgZ2VuZSBsZW5ndGhzIGZyb20gR0ZGIChzaG91bGQgYmUgaW4gYmFzZSBwYWlycykKZ2VuZV9sZW5ndGhzIDwtIGdmZl9hbm5vdCR3aWR0aCAgCm5hbWVzKGdlbmVfbGVuZ3RocykgPC0gcm93bmFtZXMoZ2ZmX2Fubm90KSAgIyBzZXQgZ2VuZSBuYW1lcwoKIyBNYWtlIHN1cmUgZ2VuZV9sZW5ndGhzIGFyZSBpbiB0aGUgc2FtZSBvcmRlciBhcyByYXdjb3VudHMgcm93cwpnZW5lX2xlbmd0aHMgPC0gZ2VuZV9sZW5ndGhzW3Jvd25hbWVzKHJhd2NvdW50cyldCgojIE1ha2UgREdFTGlzdCBvYmplY3QgZm9yIGVkZ2VSCmRnZSA8LSBER0VMaXN0KGNvdW50cyA9IHJhd2NvdW50cykKCiMgQ2FsY3VsYXRlIFJQS00gKHJlYWRzIHBlciBraWxvYmFzZSBwZXIgbWlsbGlvbikKcnBrbSA8LSBycGttKGRnZSwgZ2VuZS5sZW5ndGggPSBnZW5lX2xlbmd0aHMpCgojIFdyYXAgUlBLTSBpbiBhbiBleHBlcmltZW50IG9iamVjdCB0b28KZXhwdF9ycGttIDwtIGNyZWF0ZV9leHB0KGNvdW50X2RhdGFmcmFtZSA9IHJwa20sCiAgICAgICAgICAgICAgICAgICAgICAgICBtZXRhZGF0YSA9IG1ldGFkYXRhLCAKICAgICAgICAgICAgICAgICAgICAgICAgIGdlbmVfaW5mbyA9IGdmZl9hbm5vdCkKCiMgbG9nMiB0cmFuc2Zvcm0gdGhlIFJQS00gKGFnYWluLCBubyBub3JtYWxpemF0aW9uIGJleW9uZCBsb2cpCnJwa21fbm9ybWFsaXplZF9leHAgPC0gbm9ybWFsaXplX2V4cHQoZXhwdF9ycGttLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0cmFuc2Zvcm0gPSAibG9nMiIsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG5vcm0gPSAicmF3IikKCiMgR2V0IHRoZSBhY3R1YWwgZXhwcmVzc2lvbiBtYXRyaWNlcyBmcm9tIHRoZSBleHBlcmltZW50IG9iamVjdHMKcnBrbSA8LSBleHBycyhycGttX25vcm1hbGl6ZWRfZXhwKQpjcG0gPC0gZXhwcnMoY3BtX25vcm1hbGl6ZWRfZXhwKQoKIyBNZXJnZSBDUE0gYW5kIFJQS00gbWF0cmljZXMgd2l0aCBnZW5lIGFubm90YXRpb24gaW5mbwpjcG1fYW5ub3RhdGVkIDwtIG1lcmdlKGdmZl9hbm5vdCwgY3BtLCBieSA9IDApCnJwa21fYW5ub3RhdGVkIDwtIG1lcmdlKGdmZl9hbm5vdCwgcnBrbSwgYnkgPSAwKQpgYGAKCiMjIEhlYXRtYXAgb2YgYWxnaW5hdGUgZ2VuZXMKVGhpcyBnZW5lcmF0ZXMgYSBzbWFsbCBoZWF0bWFwIG9mIHRoZSBhbGdpbmF0ZSBnZW5lcyBhY3Jvc3Mgc2FtcGxlcyB0byB2aXN1YWxpemUgZXhwcmVzc2lvbi4KCmBgYHtyIGhlYXRtYXBzfQojIExvb3AgdGhyb3VnaCBDUE0gYW5kIFJQS00gdG8gcGxvdCBib3RoIGhlYXRtYXBzCmZvciAodHlwZSBpbiBjKCJjcG0iLCAicnBrbSIpKSB7CiAgCiAgIyBEZWZpbmUgZ2VuZSBsaXN0CiAgYWxnX29wZXJvbl9nZW5lcyA8LSBjKAogICAgImFsZ0EiLCAiYWxnRiIsICJhbGdKIiwgImFsZ0kiLCAiYWxnTCIsICJhbGdYIiwgImFsZ0ciLCAiYWxnRSIsCiAgICAiYWxnSyIsICJhbGc0NCIsICJhbGc4IiwgImFsZ0QiLCAibXVjQyIsICJtdWNCIiwgImFsZ1UiLCAiYWxnVyIsCiAgICAiYWxnUCIsICJhbGdRIiwgImFsZ1IiLCAiYWxnWiIsICJhbGdDIiwgImFsZ0IiCiAgKQoKICAjIEdldCB0aGUgcmlnaHQgZGF0YSBwZXIgdHlwZSBvZiB0cmFuc2Zvcm1hdGlvbgogIGRhdGFfYW5ub3QgPC0gZ2V0KHBhc3RlMCh0eXBlLCAiX2Fubm90YXRlZCIpKQogIAogICMgRml4IHJvd25hbWVzIHVzaW5nIGdlbmUgbmFtZSAobWFrZSB1bmlxdWUpCiAgcm93bmFtZXMoZGF0YV9hbm5vdCkgPC0gbWFrZS5uYW1lcyhkYXRhX2Fubm90JG5hbWUsIHVuaXF1ZSA9IFRSVUUpCiAgCiAgIyBEcm9wIHVuZGVzaXJlZCBtZXRhZGF0YSBjb2x1bW5zCiAgY29sc190b19leGNsdWRlIDwtIGMoIlJvdy5uYW1lcyIsInNlcW5hbWVzIiwic3RhcnQiLCJlbmQiLCJ3aWR0aCIsCiAgICAgICAgICAgICAgICAgICAgICAgInN0cmFuZCIsInNvdXJjZSIsInR5cGUiLCJzY29yZSIsInBoYXNlIiwKICAgICAgICAgICAgICAgICAgICAgICAiZ2VuZV9pZCIsIk5hbWUiLCJEYnhyZWYiLCJBbGlhcyIsIm5hbWUiLAogICAgICAgICAgICAgICAgICAgICAgICJQYXJlbnQiLCJsb2N1cyIpCiAga2VlcCA8LSAhbmFtZXMoZGF0YV9hbm5vdCkgJWluJSBjb2xzX3RvX2V4Y2x1ZGUKICBleHByX2RhdGEgPC0gZGF0YV9hbm5vdFssIGtlZXBdCiAgCiAgIyBFeHRyYWN0IG1hdHJpeCBmb3IgYWxnIGdlbmVzCiAgYWxnX2dlbmVzIDwtIGV4cHJfZGF0YVthbGdfb3Blcm9uX2dlbmVzLCAsIGRyb3AgPSBGQUxTRV0KICBtYXQgPC0gYWxnX2dlbmVzIC0gcm93TWVhbnMoYWxnX2dlbmVzKQogIAogICMgUHVsbCBzYW1wbGUgYW5ub3RhdGlvbnMKICBhbm5vIDwtIGFzLmRhdGEuZnJhbWUobWV0YWRhdGFbLCBjKCJjb25kaXRpb24iLCAiaW5kdWNlZCIpXSkKICBjb2xuYW1lcyhhbm5vKSA8LSBjKCJDb25kaXRpb24iLCAiSW5kdWNlZCIpCiAgCiAgIyBDb2xvciBzZXR0aW5ncyBmb3IgYW5ub3RhdGlvbgogIyBhbm5fY29sb3JzIDwtIGxpc3QoCiAjICBDb25kaXRpb24gPSBjKCJMQi1ETVNPIiA9ICIjRDk1RjAyIiwgIkxCLUlQVEciID0gIiM2NkE2MUUiKSwKICMgIEluZHVjZWQgPSBjKCJubyIgPSAibGlnaHRncmV5IiwgInllcyIgPSAiZ3JheTI3IikKICMgKQogIAogICMgU2NhbGUgZGF0YSBieSByb3cgKFotc2NvcmUgcGVyIGdlbmUpCiAgc2NhbGVkX21hdCA8LSB0KHNjYWxlKHQobWF0KSkpCgogICMgQnVpbGQgaGVhdG1hcAogIGh0IDwtIEhlYXRtYXAoCiAgICBzY2FsZWRfbWF0LAogICAgbmFtZSA9ICJaLXNjb3JlIiwKICAgIHRvcF9hbm5vdGF0aW9uID0gSGVhdG1hcEFubm90YXRpb24oZGYgPSBhbm5vIywgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICNjb2wgPSBhbm5fY29sb3JzCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICksCiAgICBjbHVzdGVyX3Jvd3MgPSBUUlVFLAogICAgY2x1c3Rlcl9jb2x1bW5zID0gVFJVRSwKICAgIHJvd19uYW1lc19zaWRlID0gImxlZnQiLAogICAgY29sdW1uX25hbWVzX3NpZGUgPSAiYm90dG9tIiwKICAgIGNvbHVtbl9uYW1lc19yb3QgPSA0NSwKICAgIHJvd19uYW1lc19ncCA9IGdwYXIoZm9udHNpemUgPSAxMCksCiAgICBjb2x1bW5fbmFtZXNfZ3AgPSBncGFyKGZvbnRzaXplID0gMTApLAogICAgc2hvd19yb3dfZGVuZCA9IFRSVUUsCiAgICBzaG93X2NvbHVtbl9kZW5kID0gVFJVRQogICkKICAKICAjIERyYXcgYW5kIHNhdmUKICBkcmF3KGh0KQogIAogIGNhaXJvX3BkZihwYXN0ZTAoU3lzLkRhdGUoKSwgIl9hbGdfaGVhdG1hcF8iLCB0eXBlLCAiXyIsIGV4cCwgIl8iLCB0cmVhdG1lbnQsICIudnMuIiwgY29udHJvbCwgIi5wZGYiKSkKICBkcmF3KGh0KQogIGRldi5vZmYoKQp9CgpgYGAKCiMgREVTZXEKYGBge3IgREVTZXF9CiMgQ3JlYXRlIERFU2VxRGF0YVNldApkZHMgPC0gREVTZXFEYXRhU2V0RnJvbU1hdHJpeCgKICBjb3VudERhdGEgPSByYXdjb3VudHMsCiAgY29sRGF0YSA9IG1ldGFkYXRhLAogIGRlc2lnbiA9IGRlc2lnbl9mb3JtdWxhCikKCiMgQ2hlY2sgbnVtYmVyIG9mIGdlbmVzIGJlZm9yZSBmaWx0ZXJpbmcKcHJpbnQoIk51bWJlciBvZiBnZW5lcyBiZWZvcmUgZmlsdGVyaW5nOiIpCm5yb3coZGRzKQoKIyBGaWx0ZXI6IGtlZXAgZ2VuZXMgd2l0aCBjb3VudHMgPj0xMCBpbiBhdCBsZWFzdCAyIHNhbXBsZXMKcHJpbnQoIkZpbHRlcmluZyBnZW5lcyB3aXRoIGNvdW50cyA+PTEwIGluIGF0IGxlYXN0IDIgc2FtcGxlcyIpCmtlZXAgPC0gcm93U3Vtcyhjb3VudHMoZGRzKSA+PSAxMCkgPj0gMiAKZGRzX2ZpbHRlcmVkIDwtIGRkc1trZWVwLCBdCgojIENoZWNrIG51bWJlciBvZiBnZW5lcyBhZnRlciBmaWx0ZXJpbmcKcHJpbnQoIk51bWJlciBvZiBnZW5lcyBhZnRlciBmaWx0ZXJpbmc6IikKbnJvdyhkZHNfZmlsdGVyZWQpCgojIFJ1biBERVNlcSBvbiB0aGUgZmlsdGVyZWQgZGF0YXNldApkZHNfZmluYWwgPC0gREVTZXEoZGRzX2ZpbHRlcmVkKQoKIyBTaG93IGF2YWlsYWJsZSByZXN1bHRzCnJlc3VsdHNOYW1lcyhkZHNfZmluYWwpCgojIFBsb3QgZGlzcGVyc2lvbiBlc3RpbWF0ZXMKcGxvdERpc3BFc3RzKGRkc19maW5hbCkKYGBgCgojIFNhdmUgREUgdGFibGVzIApgYGB7ciBzYXZlX3RhYmxlc30KcmVzIDwtIHJlc3VsdHMoZGRzX2ZpbmFsLCBjb250cmFzdCA9IGMoImV4cF8xMSIsICJ0cmVhdG1lbnQiLCAiY29udHJvbCIpKQpyZXNkZiA8LSBhcy5kYXRhLmZyYW1lKHJlcykKcmVzZGYkZ2VuZV9pZCA8LSByb3duYW1lcyhyZXNkZikKCmNwbSA8LSBhcy5kYXRhLmZyYW1lKGNwbSkKcnBrbSA8LSBhcy5kYXRhLmZyYW1lKHJwa20pCgojIGFwcGVuZCB0eXBlIG9mIHZhbHVlIHRvIHRoZSBkYXRhZnJhbWUgICAKY29sbmFtZXMocnBrbSkgPC0gcGFzdGUwKGNvbG5hbWVzKHJwa20pLCAiX3Jwa20iKQpjb2xuYW1lcyhjcG0pIDwtIHBhc3RlMChjb2xuYW1lcyhjcG0pLCAiX2NwbSIpCnJhd2NvdW50czIgPC0gcmF3Y291bnRzCmNvbG5hbWVzKHJhd2NvdW50czIpIDwtIHBhc3RlMChjb2xuYW1lcyhyYXdjb3VudHMyKSwgIl9yYXdjb3VudHMiKQoKIyBlbnN1cmUgY29sdW1uIGdlbmVfaWQgaXMgdGhlcmUKcnBrbSRnZW5lX2lkIDwtIHJvd25hbWVzKHJwa20pCmNwbSRnZW5lX2lkIDwtIHJvd25hbWVzKGNwbSkKcmF3Y291bnRzMiRnZW5lX2lkIDwtIHJvd25hbWVzKHJhd2NvdW50czIpCgpnZmZfYW5ub3QgPC0gZ2ZmX2Fubm90W3Jvd25hbWVzKHJlc2RmKSxdCgpkaW0ocmVzZGYpCiMgbWVyZ2Ugd2l0aCBzaWdnZW5lcwpyZXNkZiA8LSBmdWxsX2pvaW4ocmVzZGYsIGdmZl9hbm5vdCwgYnkgPSAiZ2VuZV9pZCIsIGNvcHkgPSBGQUxTRSkKaGVhZChyZXNkZikKCnBhX2dvIDwtIHJlYWQuY3N2KCIuLi8uLi9yYXctZGF0YS8yMDg5NjNfbWljcm9iZXMub25saW5lLmNzdiIsIGhlYWRlciA9IFRSVUUpCmhlYWQocGFfZ28pCgpjb2xfdG9fcmVtb3ZlIDwtIGludGVyc2VjdChjb2xuYW1lcyhwYV9nbyksY29sbmFtZXMocmVzZGYpKQpjb2xfdG9fcmVtb3ZlIDwtIHNldGRpZmYoY29sX3RvX3JlbW92ZSwgIkFsaWFzIikgIyByZXR1cm5zIGFsbCBleGNlcHQgdGhlIHNlY29uZCB2ZWN0b3IKCnBhX2dvIDwtIHBhX2dvICU+JSBzZWxlY3QoLWFsbF9vZihjb2xfdG9fcmVtb3ZlKSkKCnJlc2RmIDwtIGZ1bGxfam9pbihyZXNkZiwgcnBrbSwgYnkgPSAiZ2VuZV9pZCIsIGNvcHkgPSBGQUxTRSkKcmVzZGYgPC0gZnVsbF9qb2luKHJlc2RmLCBjcG0sIGJ5ID0gImdlbmVfaWQiLCBjb3B5ID0gRkFMU0UpCnJlc2RmIDwtIGZ1bGxfam9pbihyZXNkZiwgcmF3Y291bnRzMiwgYnkgPSAiZ2VuZV9pZCIsIGNvcHkgPSBGQUxTRSkKI3Jlc2RmIDwtIGZ1bGxfam9pbihyZXNkZiwgcGFfZ28sIGJ5ID0gIkFsaWFzIiwgY29weSA9IEZBTFNFKQpyZXNkZiA8LSBtZXJnZShyZXNkZiwgcGFfZ28sIGJ5ID0gIkFsaWFzIiwgY29weSA9IEZBTFNFKQoKaGVhZChyZXNkZikKd3JpdGUuY3N2KHJlc2RmLCBwYXN0ZTAoU3lzLkRhdGUoKSAsICJfcmVzdWx0c18iLGV4cCwgIl8iLCB0cmVhdG1lbnQsIi52cy4iLCBjb250cm9sLCIuY3N2IiksIHJvdy5uYW1lcyA9IEZBTFNFKQpgYGAKCiMgU2F2ZSBERUcgbGlzdApgYGB7ciBzYXZlX2RlZ30KIyBmaWx0ZXIgZm9yIHNpZ25pZmljYW50bHkgZGlmZmVyZW50aWFsbHkgZXhwcmVzc2VkIGdlbmVzIChMRkMgPj0gMywgcGFkaiA8IHB2YWxfY3V0b2ZmKQpkaW0ocmVzZGYpCmhlYWQocmVzZGYpCgpzdW1tYXJ5KGlzLm5hKHJlc2RmJGxvZzJGb2xkQ2hhbmdlKSkKc3VtbWFyeShpcy5uYShyZXNkZiRwYWRqKSkKCnNpZ19nZW5lcyA8LSByZXNkZlthYnMocmVzZGYkbG9nMkZvbGRDaGFuZ2UpID49IGxvZ2ZjX2N1dG9mZiAmIHJlc2RmJHBhZGogPCBwdmFsX2N1dG9mZiwgXSAgJT4lIGRyb3BfbmEobG9nMkZvbGRDaGFuZ2UpCnNpZ19nZW5lcyRkaXJlY3Rpb24gPC0gaWZlbHNlKHNpZ19nZW5lcyRsb2cyRm9sZENoYW5nZSA+IDAsICJ1cCIsICJkb3duIikKCnN1bW1hcnkoaXMubmEoc2lnX2dlbmVzJGxvZzJGb2xkQ2hhbmdlKSkKc3VtbWFyeShpcy5uYShzaWdfZ2VuZXMkcGFkaikpCgpkaW0oc2lnX2dlbmVzKQpoZWFkKHNpZ19nZW5lcykKd3JpdGUuY3N2KHNpZ19nZW5lcywgcGFzdGUwKFN5cy5EYXRlKCksIl9zaWduaWZpY2FudGdlbmVzXyIsIGV4cCwgIl8iLCB0cmVhdG1lbnQsIi52cy4iLCBjb250cm9sLCIuY3N2IiksIHJvdy5uYW1lcyA9IEZBTFNFKQpgYGAKCiMgUENBIHBsb3QgCmBgYHtyIFBDQX0KdnN0IDwtIERFU2VxRGF0YVNldEZyb21NYXRyaXgoY291bnREYXRhID0gcmF3Y291bnRzLCBjb2xEYXRhID0gbWV0YWRhdGEsIGRlc2lnbiA9IH4gZXhwXzExICsgYmF0Y2gpCgojIFZhcmlhbmNlIHN0YWJpbGl6aW5nIHRyYW5zZm9ybWF0aW9ucwp2c2QgPC0gdnN0KHZzdCwgYmxpbmQgPSBUUlVFKQoKIyBQQ0EgcGxvdApwY2FEYXRhIDwtIHBsb3RQQ0EodnNkLCBpbnRncm91cCA9IGMoImV4cF8xMSIsICJiYXRjaCIsICJjb25kaXRpb24iKSwgcmV0dXJuRGF0YSA9IFRSVUUpCnBlcmNlbnRWYXIgPC0gcm91bmQoMTAwICogYXR0cihwY2FEYXRhLCAicGVyY2VudFZhciIpKQoKIyBDYWxjdWxhdGUgdGhlIG1heGltdW0gYWJzb2x1dGUgdmFsdWUgYW1vbmcgUEMxIGFuZCBQQzIKbWF4X2Fic192YWwgPC0gbWF4KGMoYWJzKG1pbihwY2FEYXRhJFBDMSkpLCBhYnMobWluKHBjYURhdGEkUEMyKSksIGFicyhtYXgocGNhRGF0YSRQQzEpKSwgYWJzKG1heChwY2FEYXRhJFBDMikpKSkKCmF4aXNfc2l6ZSA8LSAxMAp0aGVtZV9zZXQodGhlbWVfbWluaW1hbChiYXNlX2ZhbWlseSA9ICJBcmlhbCIpKQoKcGNhIDwtIGdncGxvdChwY2FEYXRhLCBhZXMoUEMxLCBQQzIsIGNvbG9yID0gY29uZGl0aW9uLCBzaGFwZSA9IGNvbmRpdGlvbikpICsKICBnZW9tX3BvaW50KHNpemUgPSAzKSArCiAgeGxhYihwYXN0ZTAoIlBDMTogIiwgcGVyY2VudFZhclsxXSwgIiUgdmFyaWFuY2UiKSkgKwogIHlsYWIocGFzdGUwKCJQQzI6ICIsIHBlcmNlbnRWYXJbMl0sICIlIHZhcmlhbmNlIikpICsKICBzY2FsZV9jb2xvdXJfYnJld2VyKHBhbGV0dGUgPSAiRGFyazIiKSArCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcyA9IGMoMSwgMiwgMywgNCwgNSwgNiwgNywgOCwgOSwgMTAsIDExLCAxMiwgMTMsIDE0LCAxNSwgMTYpKSArCiAgY29vcmRfZml4ZWQocmF0aW8gPSAxLCB4bGltID0gYygtbWF4X2Fic192YWwsIG1heF9hYnNfdmFsKSwgeWxpbSA9IGMoLW1heF9hYnNfdmFsLCBtYXhfYWJzX3ZhbCkpICsKICBsYWJzKGNvbG9yID0gImNvbmRpdGlvbiIsIHNoYXBlID0gImNvbmRpdGlvbiIpICsKICBndWlkZXMoc2hhcGUgPSBndWlkZV9sZWdlbmQodGl0bGUgPSAiQ29uZGl0aW9uIiksICAjIEFkanVzdCBzaXplIG9mIHNoYXBlIGxlZ2VuZAogICAgY29sb3IgPSBndWlkZV9sZWdlbmQodGl0bGUgPSAiQ29uZGl0aW9uIikpICsKICB0aGVtZSgKICAgIHBhbmVsLmdyaWQgPSBlbGVtZW50X2JsYW5rKCksCiAgICBwYW5lbC5iYWNrZ3JvdW5kID0gZWxlbWVudF9ibGFuaygpLAogICAgYXhpcy50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSBheGlzX3NpemUpLAogICAgYXhpcy50aWNrcyA9IGVsZW1lbnRfbGluZShsaW5ld2lkdGggPSAwLjIpLAogICAgIGF4aXMubGluZSA9IGVsZW1lbnRfbGluZShsaW5ld2lkdGggPSAwLCBjb2xvciA9ICJibGFjayIpLCAgIyBTZXQgdGhlIGF4aXMgbGluZSBsaW5ld2lkdGgKICAgcGFuZWwuYm9yZGVyID0gZWxlbWVudF9yZWN0KGZpbGwgPSBOQSwgbGluZXdpZHRoID0gMC41KSwgICMgU2V0IHRoZSBwYW5lbCBib3JkZXIgbGluZXdpZHRoCiAgbGVnZW5kLmtleSA9IGVsZW1lbnRfcmVjdChmaWxsID0gIndoaXRlIiwgY29sb3IgPSBOQSkgICMgU2V0IGxlZ2VuZCBrZXkgYmFja2dyb3VuZCB0byB3aGl0ZSBhbmQgcmVtb3ZlIGJvcmRlcgogICkKCnByaW50KHBjYSkKCmNhaXJvX3BkZihwYXN0ZTAoU3lzLkRhdGUoKSwiX3BjYV8iLGV4cCwgIl8iLCB0cmVhdG1lbnQsIi52cy4iLCBjb250cm9sLCIucGRmIiwgc2VwID0gIiIpKQpwY2EKZGV2Lm9mZigpCgpgYGAKCiMgVm9sY2FubyBwbG90CmBgYHtyIHZvbGNhbm99CnNlcXVlbmNlIDwtIGMoInRyZWF0bWVudCIpCgpmb3IgKGkgaW4gc2VxdWVuY2UpIHsKICAgIHNldC5zZWVkKDApCiAgICByZXMgPC0gcmVzdWx0cyhkZHNfZmluYWwsIGNvbnRyYXN0ID0gYygiZXhwXzExIiwgaSwgImNvbnRyb2wiKSkKCiAgZ3AyIDwtIHdlc19wYWxldHRlcyRaaXNzb3UxCiAga2V5IDwtIHJlYWQuY3N2KCIuLi8uLi9yYXctZGF0YS9wYTE0X2dlbmVfa2V5LmNzdiIsICBoZWFkZXIgPSBUUlVFKQoKICByb3duYW1lcyhrZXkpIDwtIGtleSRnZW5lX2lkCgogIHJlcyA8LSBtZXJnZShrZXksIGFzLmRhdGEuZnJhbWUocmVzKSwgYnk9MCkKICByZXMgPC0gYXMuZGF0YS5mcmFtZShyZXMpCiAgcmVzJEFsaWFzIDwtIHJvd25hbWVzKHJlcykKCiAjIERlZmluZSBnZW5lIGxpc3QKICBhbGdfZ2VuZXMgPC0gYygKICAgICJhbGdBIiwgImFsZ0YiLCAiYWxnSiIsICJhbGdJIiwgImFsZ0wiLCAiYWxnWCIsICJhbGdHIiwgImFsZ0UiLAogICAgImFsZ0siLCAiYWxnNDQiLCAiYWxnOCIsICJhbGdEIiwgIm11Y0MiLCAibXVjQiIsICJhbGdVIiwgImFsZ1ciLAogICAgImFsZ1AiLCAiYWxnUSIsICJhbGdSIiwgImFsZ1oiLCAiYWxnQyIsICJhbGdCIgogICkKCmxvZ2ZjX2N1dG9mZgprZXl2YWxzIDwtIGlmZWxzZSgKIHJlcyRsb2cyRm9sZENoYW5nZSA8PSAtbG9nZmNfY3V0b2ZmICYgcmVzJHBhZGogPD0gcHZhbF9jdXRvZmYsIGdwMlsxXSwKICBpZmVsc2UocmVzJGxvZzJGb2xkQ2hhbmdlID49IGxvZ2ZjX2N1dG9mZiAmIHJlcyRwYWRqIDw9IHB2YWxfY3V0b2ZmLCBncDJbNV0sCiAgICAgICAgICdsaWdodGdyZXknKSkKa2V5dmFsc1tpcy5uYShrZXl2YWxzKV0gPC0gJ2xpZ2h0Z3JleScKbmFtZXMoa2V5dmFscylba2V5dmFscyA9PSBncDJbNV1dIDwtICdVcC1yZWd1bGF0ZWQnCm5hbWVzKGtleXZhbHMpW2tleXZhbHMgPT0gJ2xpZ2h0Z3JleSddIDwtICdOUycKbmFtZXMoa2V5dmFscylba2V5dmFscyA9PSBncDJbMV1dIDwtICdEb3duLXJlZ3VsYXRlZCcKCnZvbGNhbm8gPC0gRW5oYW5jZWRWb2xjYW5vKHJlcywKICAgICAgICAgICAgICAgICAgICAgbGFiID0gcmVzW1sibmFtZSJdXSwKICAgICAgICAgICAgICAgICAgICAgc2VsZWN0TGFiID0gYWxnX2dlbmVzLAogICAgICAgICAgICAgICAgICAgICB4ID0gJ2xvZzJGb2xkQ2hhbmdlJywKICAgICAgICAgICAgICAgICAgICAgeSA9ICdwYWRqJywKICAgICAgICAgICAgICAgICAgICAgeGxhYiA9IGJxdW90ZSh+TG9nWzJdfiAnZm9sZCBjaGFuZ2UnKSwKICAgICAgICAgICAgICAgICAgICAgcEN1dG9mZiA9IHB2YWxfY3V0b2ZmLAogICAgICAgICAgICAgICAgICAgICBGQ2N1dG9mZiA9IGxvZ2ZjX2N1dG9mZiwKICAgICAgICAgICAgICAgICAgICAgY3V0b2ZmTGluZVR5cGUgPSAiZGFzaGVkIiwKICAgICAgICAgICAgICAgICAgICAgcG9pbnRTaXplID0gMSwKICAgICAgICAgICAgICAgICAgICAgbGFiU2l6ZSA9IDMuMCwKICAgICAgICAgICAgICAgICAgICAgYXhpc0xhYlNpemUgPSAxNCwKICAgICAgICAgICAgICAgICAgICAgbGFiQ29sID0gJ2JsYWNrJywKICAgICAgICAgICAgICAgICAgICAgbGFiRmFjZSA9ICdib2xkJywKICAgICAgICAgICAgICAgICAgICAgI2NvbCA9IGMoJ2dyZXknLCBncDJbMl0sIGdwMlszXSwgZ3AyWzFdKSwKICAgICAgICAgICAgICAgICAgICAgY29sQ3VzdG9tID0ga2V5dmFscywKICAgICAgICAgICAgICAgICAgICAgY29sQWxwaGEgPSA1LzUsCiAgICAgICAgICAgICAgICAgICAgIGxlZ2VuZFBvc2l0aW9uID0gJ2JvdHRvbScsCiAgICAgICAgICAgICAgICAgICAgIGxlZ2VuZExhYlNpemUgPSAxMiwKICAgICAgICAgICAgICAgICAgICAgbGVnZW5kSWNvblNpemUgPSA0LjAsCiAgICAgICAgICAgICAgICAgICAgIGRyYXdDb25uZWN0b3JzID0gVFJVRSwKICAgICAgICAgICAgICAgICAgICAgYXJyb3doZWFkcyA9IEZBTFNFLAogICAgICAgICAgICAgICAgICAgICB3aWR0aENvbm5lY3RvcnMgPSAwLjIsCiAgICAgICAgICAgICAgICAgICAgIGdyaWRsaW5lcy5tYWpvciA9IEZBTFNFLAogICAgICAgICAgICAgICAgICAgICBib3JkZXJXaWR0aCA9IDAuMywKICAgICAgICAgICAgICAgICAgICAgZ3JpZGxpbmVzLm1pbm9yID0gRkFMU0UsCiAgICAgICAgICAgICAgICAgICAgIHRpdGxlID0gcGFzdGUwKHRyZWF0bWVudCwiIHZzICIsIGNvbnRyb2wsIHNlcCA9ICIiKSwKICAgICAgICAgICAgICAgICAgICAgc3VidGl0bGUgPSAiIiwKICAgICAgICAgICAgICAgICAgICAgI2NhcHRpb24gPSBwYXN0ZTAoInRvdGFsID0gIiwgbnJvdyh0b3B0YWJsZSksICIgdmFyaWFibGVzIiksCiAgICAgICAgICAgICAgICAgICAgIGNhcHRpb24gPSBwYXN0ZTAoIlZvbGNhbm8gcGxvdCB3aXRoIGxvZzIgRkMgb2YgIiwgbG9nZmNfY3V0b2ZmLCIgYW5kIFAgdmFsdWUgb2YgIixwdmFsX2N1dG9mZiwgIi4iKSwKICAgICAgICAgICAgICAgICAgICAgbWF4Lm92ZXJsYXBzID0gMzAsCiAgICAgICAgICAgICAgICAgICAgICN5bGltID0gYygwLDEwKSwKICAgICAgICAgICAgICAgICAgICAgY29sQ29ubmVjdG9ycyA9ICdibGFjaycpICsgCiAgICAgIHRoZW1lKGF4aXMudGlja3M9ZWxlbWVudF9saW5lKGxpbmV3aWR0aCA9MC4yKSkgIyArCiAjIHNjYWxlX3hfY29udGludW91cyhicmVha3MgPSBjKC01LCAtMi41LCAwLCAyLjUsIDUpLCBsaW1pdHMgPSBjKC02LDYpKQoKcHJpbnQodm9sY2FubykKCmNhaXJvX3BkZihwYXN0ZTAoU3lzLkRhdGUoKSwiX3ZvbGNhbm9fIiwgZXhwLCAiXyIsIHRyZWF0bWVudCwiLnZzLiIsIGNvbnRyb2wsIi5wZGYiLCBzZXAgPSAiIikpCnByaW50KHZvbGNhbm8pCmRldi5vZmYoKQp9CgpgYGAKCiMgTUEgcGxvdApgYGB7ciBNQV9wbG90fQpzZXF1ZW5jZSA8LSBjKCJ0cmVhdG1lbnQiKQoKZm9yIChpIGluIHNlcXVlbmNlKSB7CiAgICBzZXQuc2VlZCgwKQoKcmVzIDwtIHJlc3VsdHMoZGRzX2ZpbmFsLCBjb250cmFzdCA9IGMoImV4cF8xMSIsIGksICJjb250cm9sIikpCgpncDIgPC0gd2VzX3BhbGV0dGVzJFppc3NvdTEKCnJlcyRnZW5lX2lkIDwtIHJvd25hbWVzKHJlcykKcmVzIDwtIGFzLmRhdGEuZnJhbWUocmVzKQpoZWFkKHJlcykKcmVzIDwtIGZ1bGxfam9pbihrZXksIGFzLmRhdGEuZnJhbWUocmVzKSwgYnk9ImdlbmVfaWQiKQpyb3duYW1lcyhyZXMpIDwtIHJlcyRBbGlhcwoKYW55KGlzLm5hKHJlcykpICAjIENoZWNrIGlmIHRoZXJlIGFyZSBhbnkgTkFzCgpyZXMgPC0gcmVzICU+JSBmaWx0ZXIoIWlzLm5hKGJhc2VNZWFuKSkKCgptYSA8LSBnZ21hcGxvdChyZXMsIAogICAgICAgICAgICAgICBmZHIgPSBwdmFsX2N1dG9mZiwgCiAgICAgICAgICAgICAgIGZjID0gMl4obG9nZmNfY3V0b2ZmKSwgCiAgICAgICAgICAgICAgIHNpemUgPSAxLAogICAgICAgICAgICAgICBwYWxldHRlID0gYygoZ3AyWzVdKSwgKGdwMlsxXSkpLAogICAgICAgICAgICAgICBzZWxlY3QudG9wLm1ldGhvZCA9ICJmYyIsCiAgICAgICAgICAgICAjICB5bGltID0gYygtMTUsMTUpLAogICAgICAgICAgICAgICBsZWdlbmQgPSAiYm90dG9tIiwgCiAgICAgICAgICAgICAgIHRpdGxlID0gcGFzdGUwKHRyZWF0bWVudCwiIHZzICIsIGNvbnRyb2wsIHNlcCA9ICIiKSwKICAgICAgICAgICAgICBzdWJ0aXRsZSA9IHBhc3RlMCgiVm9sY2FubyBwbG90IHdpdGggbG9nMiBGQyBvZiAiLCBsb2dmY19jdXRvZmYsIiBhbmQgUCB2YWx1ZSBvZiAiLHB2YWxfY3V0b2ZmLCAiLiIpLAogICAgICAgICAgICAgICBmb250LmxhYmVsID0gYyggMTIpLCAKICAgICAgICAgICAgIGxhYmVsLnJlY3RhbmdsZSA9IEZBTFNFLAogICAgICAgICAgICAgIGdlbmVuYW1lcyA9IHJlcyRuYW1lLAogICAgICAgICAgICAgIGZvbnQubGVnZW5kID0gYygxMiksCiAgICAgICAgICAgICAgIGZvbnQubWFpbiA9IGMoImJvbGQiLCAxNCksCiAgICAgICAgICAgICAgIGFscGhhID0gMS41LAogICAgICAgICAgICAgICBnZ3RoZW1lID0gZ2dwbG90Mjo6dGhlbWUocGFuZWwuZ3JpZCA9ZWxlbWVudF9ibGFuaygpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgYXhpcy50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSAxMiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBheGlzLnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSAxNCksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBwYW5lbC5iYWNrZ3JvdW5kID0gZWxlbWVudF9ibGFuaygpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgcGFuZWwuYm9yZGVyID0gZWxlbWVudF9yZWN0KGZpbGwgPSBOQSkpKSAKcHJpbnQobWEpCgpjYWlyb19wZGYocGFzdGUwKFN5cy5EYXRlKCksIl9tYV8iLCBleHAsICJfIiwgdHJlYXRtZW50LCIudnMuIiwgY29udHJvbCwiLnBkZiIpKQpwcmludChtYSkKZGV2Lm9mZigpCn0KYGBgCgo=