Human Macrophages: atb version of M-CSF v GM-CSF: LPS, LPS+Adenosine, LPS+PGE2

In this version of the worksheet, I am hoping to perform basically the same analyses, but do it in a fashion which is more in my own style.

1 Annotation data

Collect the human annotation data using biomaRt.

gene_info <- load_biomart_annotations(host="useast.ensembl.org")$annotation

## The biomart annotations file already exists, loading from it.

rownames(gene_info) <- make.names(gene_info[["ensembl_gene_id"]], unique=TRUE)
tx_gene_map <- gene_info[, c("ensembl_transcript_id", "ensembl_gene_id")]

2 Experimental design

I am going to use Kajal’s sample sheet without modification.

design <- read.table("sample_sheets/MetaData only 4 hour.txt", header=TRUE, sep='\t')
design[["Patient"]] <- as.factor(design[["Patient"]])
design[["Stimulation"]] <- as.factor(design[["Stimulation"]])
design[["Batch"]] <- as.factor(design[["Batch"]])
design[["Growth"]] <- as.factor(design[["Growth"]])
files <- file.path("kallisto abundance files/", design$HPGL.Identifier, "abundance.tsv")
names(files) <- paste0("HPGL09", c(12:31, 42:60))
rownames(design) <- design[[1]]
design[["condition"]] <- design[["Stimulation"]]
design[["file"]] <- glue::glue("preprocessing/{rownames(design)}/abundance.tsv")
colnames(design) <- tolower(colnames(design))
design[["gp"]] <- as.factor(glue("{design[['growth']]}_{design[['stimulation']]}"))

## Set up a column called stim_pred which is a predicate of stimulated vs. unstimulated samples.
design[["stim_pred"]] <- "stimulated"
ns_idx <- design[["stimulation"]] == "NS"
design[ns_idx, "stim_pred"] <- "unstimulated"

3 Create expressionset

We have some annotation data and experimental metadata.

hs_expt <- create_expt(metadata=design, gene_info=gene_info, tx_gene_map=tx_gene_map)

## Reading the sample metadata.

## The sample definitions comprises: 39 rows(samples) and 10 columns(metadata fields).

## Reading count tables.

## Using the transcript to gene mapping.

## Reading kallisto data with tximport.

## Finished reading count data.

## Matched 34431 annotations and counts.

## Bringing together the count matrix and gene information.

## The mapped IDs are not the rownames of your gene information, changing them now.

## Some annotations were lost in merging, setting them to 'undefined'.

## Saving the expressionset to 'expt.rda'.

## The final expressionset has 34431 rows and 39 columns.

hs_expt <- set_expt_batches(hs_expt, fact="growth")
stim_expt <- subset_expt(hs_expt, subset="stimulation!='NS'")

## Using a subset expression.

## There were 39, now there are 30 samples.

4 Write out the expressionset data

We can write out the data to an excel file in the hopes that it will prove useful.

written <- write_expt(hs_expt, batch="raw",
                      excel=glue::glue("excel/hs_expt-v{ver}.xlsx"))

## Writing the first sheet, containing a legend and some summary data.

## Writing the raw reads.

## Graphing the raw reads.

## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'

## Attempting mixed linear model with: ~  (1|condition) + (1|batch)

## Fitting the expressionset to the model, this is slow.

## Dividing work into 100 chunks...

## 
## Total:163 s

## Placing factor: condition at the beginning of the model.

## Writing the normalized reads.

## Graphing the normalized reads.

## `geom_smooth()` using formula 'y ~ x'

## `geom_smooth()` using formula 'y ~ x'

## Attempting mixed linear model with: ~  (1|condition) + (1|batch)

## Fitting the expressionset to the model, this is slow.

## Dividing work into 100 chunks...

## 
## Total:154 s

## Placing factor: condition at the beginning of the model.

## Writing the median reads by factor.

There is an important caveat, this is not taking into account the patient effects.

4.1 Show some written plots

written$legend_plot

written$raw_libsize

written$raw_density

written$raw_boxplot

written$norm_nonzero

written$norm_corheat

written$norm_disheat

written$norm_pca

5 Consider different models for the data

We may wish to lower the variance from the patients and/or the GM/M effects.

5.1 Current state with sva

hs_batch <- normalize_expt(hs_expt, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="svaseq")

## This function will replace the expt$expressionset slot with:

## log2(svaseq(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Warning in normalize_expt(hs_expt, transform = "log2", convert = "cpm", :
## Quantile normalization and sva do not always play well together.

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with svaseq.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 6 surrogate variable(s).

## Attempting svaseq estimation with 6 surrogates.

## There are 466 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

5.2 Current state with residual-based batch adjustment

Because we are explicitly removing the effect of GM/M, the patient effect really becomes apparent.

hs_batch <- normalize_expt(hs_expt, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="limmaresid")

## This function will replace the expt$expressionset slot with:

## log2(limmaresid(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with limmaresid.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 6 surrogate variable(s).

## batch_counts: Using residuals of limma's lmfit to remove batch effect.

## There are 236935 (48%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

5.3 Set batch to patient and repeat sva

The picture with sva should be the same as the first plot, just with 6 shapes instead of two.

hs_pat <- set_expt_batches(hs_expt, fact="patient")
hs_batch <- normalize_expt(hs_pat, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="svaseq")

## This function will replace the expt$expressionset slot with:

## log2(svaseq(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Warning in normalize_expt(hs_pat, transform = "log2", convert = "cpm", norm =
## "quant", : Quantile normalization and sva do not always play well together.

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with svaseq.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 6 surrogate variable(s).

## Attempting svaseq estimation with 6 surrogates.

## There are 466 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

5.4 Patient batch and limma

This picture should be different, and should show us the M/GM effect as opposed to the patient effect.

hs_batch <- normalize_expt(hs_pat, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="limmaresid")

## This function will replace the expt$expressionset slot with:

## log2(limmaresid(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with limmaresid.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 6 surrogate variable(s).

## batch_counts: Using residuals of limma's lmfit to remove batch effect.

## There are 239301 (49%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

5.5 Combine growth and stimulation and repeat sva

I am not sure what this will look like. Since two aspects of the data are in the condition portion of the model matrix, I think it should look different. When I created the design matrix, I made a column for this purpose; I called it ‘gp’, but honestly I don’t remember why…

hs_gs <- set_expt_conditions(hs_pat, fact="gp")
hs_batch <- normalize_expt(hs_gs, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="svaseq")

## This function will replace the expt$expressionset slot with:

## log2(svaseq(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Warning in normalize_expt(hs_gs, transform = "log2", convert = "cpm", norm =
## "quant", : Quantile normalization and sva do not always play well together.

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with svaseq.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 5 surrogate variable(s).

## Attempting svaseq estimation with 5 surrogates.

## There are 640 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

It seems to me that is primarily showing us differences between M/GM.

What about limma?

hs_batch <- normalize_expt(hs_gs, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="limma")

## This function will replace the expt$expressionset slot with:

## log2(limma(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 841 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with limma.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 459917 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 841 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 28536 entries are 0<x<1: 6%.

## The be method chose 5 surrogate variable(s).

## batch_counts: Using limma's removeBatchEffect to remove batch effect.

## If you receive a warning: 'NANs produced', one potential reason is that the data was quantile normalized.

## There are 762 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

Same deal, just more. This might be the moment to reconsider the fact that Kajal’s work focused only on the M/GM samples and AFAICT ignored the non-stimulated samples.

tmp <- subset_expt(hs_gs, subset='stimulation!="NS"')

## Using a subset expression.

## There were 39, now there are 30 samples.

hs_batch <- normalize_expt(tmp, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="sva")

## This function will replace the expt$expressionset slot with:

## log2(sva(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Warning in normalize_expt(tmp, transform = "log2", convert = "cpm", norm =
## "quant", : Quantile normalization and sva do not always play well together.

## Step 1: performing count filter with option: cbcb

## Removing 22066 low-count genes (12365 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 644 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with sva.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 351625 entries are x>1: 95%.

## batch_counts: Before batch/surrogate estimation, 644 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 18681 entries are 0<x<1: 5%.

## The be method chose 5 surrogate variable(s).

## Estimate type 'sva' is shorthand for 'sva_unsupervised'.

## Other sva options include: sva_supervised and svaseq.

## Attempting sva unsupervised surrogate estimation with 5 surrogates.

## There are 338 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

tmp <- subset_expt(hs_gs, subset='stimulation!="NS"')

## Using a subset expression.

## There were 39, now there are 30 samples.

hs_batch <- normalize_expt(tmp, transform="log2", convert="cpm",
                           norm="quant", filter=TRUE, batch="ruv_empirical")

## This function will replace the expt$expressionset slot with:

## log2(ruv_empirical(cpm(quant(cbcb(data)))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Step 1: performing count filter with option: cbcb

## Removing 22066 low-count genes (12365 remaining).

## Step 2: normalizing the data with quant.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 644 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with ruv_empirical.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 351625 entries are x>1: 95%.

## batch_counts: Before batch/surrogate estimation, 644 entries are x==0: 0%.

## batch_counts: Before batch/surrogate estimation, 18681 entries are 0<x<1: 5%.

## The be method chose 5 surrogate variable(s).

## Attempting ruvseq empirical surrogate estimation with 5 surrogates.

## There are 219 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

plot_pca(hs_batch)$plot

Interesting, I did not put them all into the test block above, but I tried out a bunch of small changes to the model and adjusters. I think I learned one primary lesson: patient number 3 is a bit weird. I think I might suggest removing this person from the data.

The next lesson I learned is that LPS is way different than LA/LP, something which I kind of knew from other work, but worth remembering.

6 Differential expression

There are a few ways to consider differential expression for this data. In all cases I think it is safe to assume that we wish to use patient as the batch factor/surrogate variable.

With that in mind, here are the factors of the data to which we have usable variance/experimental design:

1. Stimulated vs. Unstimulated: This I think has the most variance in the data, even including patient. We can access this by putting stimulation state (yes/no) in the model and just running sva against everything else, or by having a model like “~ stimulation_binary + patient + growth”
2. Stimulation types: If we want to consider all LPS vs. all LP etc… we can do that in a similar fashion, by either putting stimulation state (LPS/LP/etc) in the model and just running sva. Conversely we could do “~ stimulation_state + patient + growth” in the model.
3. Growth condition: Ibid, except “~ growth + patient + stimulation”
4. Growth+Stimulation: This is the focus of Kajal’s worksheet I think and may be repeated with “~ gp + patient”

Before I run these, lets look at the variance in the data and make sure I am not full of crap.

hs_vpin <- normalize_expt(hs_expt, convert="cpm", filter=TRUE)

## This function will replace the expt$expressionset slot with:

## cpm(cbcb(data))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Leaving the data in its current base format, keep in mind that
##  some metrics are easier to see when the data is log2 transformed, but
##  EdgeR/DESeq do not accept transformed data.

## Leaving the data unnormalized.  This is necessary for DESeq, but
##  EdgeR/limma might benefit from normalization.  Good choices include quantile,
##  size-factor, tmm, etc.

## Not correcting the count-data for batch effects.  If batch is
##  included in EdgerR/limma's model, then this is probably wise; but in extreme
##  batch effects this is a good parameter to play with.

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: not normalizing the data.

## Step 3: converting the data with cpm.

## Step 4: not transforming the data.

## Step 5: not doing batch correction.

hs_varpart <- simple_varpart(hs_vpin, factors=c("stimulation", "growth", "patient"),
                             chosen_factor="patient", do_fit=TRUE)

## Attempting mixed linear model with: ~  (1|stimulation) + (1|growth) + (1|patient)

## Fitting the expressionset to the model, this is slow.

## Dividing work into 100 chunks...

## 
## Total:132 s

## Placing factor: stimulation at the beginning of the model.

## Memory usage to store result: >450.6 Mb

## Dividing work into 100 chunks...

## 
## Total:341 s

hs_varpart$partition_plot

top_40_stimulation <- hs_varpart$percent_plot
## Now show a variance boxplot for the chosen batch factor (patient)
hs_varpart$stratify_batch_plot

hs_varpart$stratify_condition_plot

percent_growth <- replot_varpart_percent(hs_varpart, column="growth")
percent_growth$plot

percent_patient <- replot_varpart_percent(hs_varpart, column="patient")
percent_patient$plot

percent_unknown <- replot_varpart_percent(hs_varpart, column="Residuals")
percent_unknown$plot

## Perform some DE

6.0.1 Patient as batch, growth and stimulation as condition, sva

hs_filt <- normalize_expt(hs_expt, filter=TRUE)

## This function will replace the expt$expressionset slot with:

## cbcb(data)

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Leaving the data in its current base format, keep in mind that
##  some metrics are easier to see when the data is log2 transformed, but
##  EdgeR/DESeq do not accept transformed data.

## Leaving the data unconverted.  It is often advisable to cpm/rpkm
##  the data to normalize for sampling differences, keep in mind though that rpkm
##  has some annoying biases, and voom() by default does a cpm (though hpgl_voom()
##  will try to detect this).

## Leaving the data unnormalized.  This is necessary for DESeq, but
##  EdgeR/limma might benefit from normalization.  Good choices include quantile,
##  size-factor, tmm, etc.

## Not correcting the count-data for batch effects.  If batch is
##  included in EdgerR/limma's model, then this is probably wise; but in extreme
##  batch effects this is a good parameter to play with.

## Step 1: performing count filter with option: cbcb

## Removing 21885 low-count genes (12546 remaining).

## Step 2: not normalizing the data.

## Step 3: not converting the data.

## Step 4: not transforming the data.

## Step 5: not doing batch correction.

pat_gs_sva <- set_expt_conditions(hs_filt, fact="gp")
pat_gs_sva <- set_expt_batches(pat_gs_sva, fact="patient")
pat_gs_sva_de <- all_pairwise(pat_gs_sva, model_batch="sva")

## batch_counts: Before batch/surrogate estimation, 485493 entries are x>1: 99%.

## batch_counts: Before batch/surrogate estimation, 2575 entries are x==0: 1%.

## batch_counts: Before batch/surrogate estimation, 413 entries are 0<x<1: 0%.

## The be method chose 6 surrogate variable(s).

## Estimate type 'sva' is shorthand for 'sva_unsupervised'.

## Other sva options include: sva_supervised and svaseq.

## Attempting sva unsupervised surrogate estimation with 6 surrogates.

## Plotting a PCA before surrogates/batch inclusion.

## Not putting labels on the plot.

## Using sva to visualize before/after batch inclusion.

## Performing a test normalization with: raw

## This function will replace the expt$expressionset slot with:

## log2(sva(cpm(cbcb(data))))

## It will save copies of each step along the way
##  in expt$normalized with the corresponding libsizes. Keep libsizes in mind
##  when invoking limma.  The appropriate libsize is non-log(cpm(normalized)).
##  This is most likely kept at:
##  'new_expt$normalized$intermediate_counts$normalization$libsizes'
##  A copy of this may also be found at:
##  new_expt$best_libsize

## Leaving the data unnormalized.  This is necessary for DESeq, but
##  EdgeR/limma might benefit from normalization.  Good choices include quantile,
##  size-factor, tmm, etc.

## Step 1: performing count filter with option: cbcb

## Removing 0 low-count genes (12546 remaining).

## Step 2: not normalizing the data.

## Step 3: converting the data with cpm.

## Step 4: transforming the data with log2.

## transform_counts: Found 2575 values equal to 0, adding 1 to the matrix.

## Step 5: doing batch correction with sva.

## Note to self:  If you get an error like 'x contains missing values' The data has too many 0's and needs a stronger low-count filter applied.

## Passing off to all_adjusters.

## batch_counts: Before batch/surrogate estimation, 458520 entries are x>1: 94%.

## batch_counts: Before batch/surrogate estimation, 2575 entries are x==0: 1%.

## batch_counts: Before batch/surrogate estimation, 28199 entries are 0<x<1: 6%.

## The be method chose 5 surrogate variable(s).

## Estimate type 'sva' is shorthand for 'sva_unsupervised'.

## Other sva options include: sva_supervised and svaseq.

## Attempting sva unsupervised surrogate estimation with 5 surrogates.

## There are 765 (0%) elements which are < 0 after batch correction.

## Setting low elements to zero.

## Not putting labels on the plot.

## Finished running DE analyses, collecting outputs.

## Comparing analyses.

keepers <- list(
    ## GM against unstimulated
    "GM_LPS_vs_GM_NS" = c("GM_LPS", "GM_NS"),
    "GM_LP_vs_GM_NS" = c("GM_LP", "GM_NS"),
    "GM_LA_vs_GM_NS" = c("GM_LA", "GM_NS"),
    ## M against unstimulated
    "M_LPS_vs_M_NS" = c("M_LPS", "M_NS"),
    "M_LP_vs_M_NS" = c("M_LP", "M_NS"),
    "M_LA_vs_M_NS" = c("M_LA", "M_NS"),
    ## GM against LPS
    "GM_LP_vs_GM_LPS" = c("GM_LP", "GM_LPS"),
    "GM_LA_vs_GM_LPS" = c("GM_LA", "GM_LPS"),
    ## M against LPS
    "M_LP_vs_M_LPS" = c("M_LP", "M_LPS"),
    "M_LA_vs_M_LPS" = c("M_LA", "M_LPS"),
    ## GM, LA vs LP
    "GM_LP_vs_GM_LA" = c("GM_LP", "GM_LA"),
    ## M, LA vs LP
    "M_LP_vs_M_LA" = c("M_LP", "M_LA"),
    ## Last, each M vs GM
    "GM_NS_vs_M_NS" = c("GM_NS", "M_NS"),
    "GM_LPS_vs_M_LPS" = c("GM_LPS", "M_LPS"),
    "GM_LP_vs_M_LP" = c("GM_LP", "M_LP"),
    "GM_LA_vs_M_LA" = c("GM_LA", "M_LA"))
pat_gs_sva_tables <- combine_de_tables(
    pat_gs_sva_de, keepers=keepers,
    excel=glue::glue("excel/pat_gs_sva_tables-v{ver}.xlsx"))

## Deleting the file excel/pat_gs_sva_tables-v20200330.xlsx before writing the tables.

## Writing a legend of columns.

## Printing a pca plot before/after surrogates/batch estimation.

## Working on 1/16: GM_LPS_vs_GM_NS which is: GM_LPS/GM_NS.

## Found inverse table with GM_NS_vs_GM_LPS

## The ebseq table is null.

## Working on 2/16: GM_LP_vs_GM_NS which is: GM_LP/GM_NS.

## Found inverse table with GM_NS_vs_GM_LP

## The ebseq table is null.

## Working on 3/16: GM_LA_vs_GM_NS which is: GM_LA/GM_NS.

## Found inverse table with GM_NS_vs_GM_LA

## The ebseq table is null.

## Working on 4/16: M_LPS_vs_M_NS which is: M_LPS/M_NS.

## Found inverse table with M_NS_vs_M_LPS

## The ebseq table is null.

## Working on 5/16: M_LP_vs_M_NS which is: M_LP/M_NS.

## Found inverse table with M_NS_vs_M_LP

## The ebseq table is null.

## Working on 6/16: M_LA_vs_M_NS which is: M_LA/M_NS.

## Found inverse table with M_NS_vs_M_LA

## The ebseq table is null.

## Working on 7/16: GM_LP_vs_GM_LPS which is: GM_LP/GM_LPS.

## Found inverse table with GM_LPS_vs_GM_LP

## The ebseq table is null.

## Working on 8/16: GM_LA_vs_GM_LPS which is: GM_LA/GM_LPS.

## Found inverse table with GM_LPS_vs_GM_LA

## The ebseq table is null.

## Working on 9/16: M_LP_vs_M_LPS which is: M_LP/M_LPS.

## Found inverse table with M_LPS_vs_M_LP

## The ebseq table is null.

## Working on 10/16: M_LA_vs_M_LPS which is: M_LA/M_LPS.

## Found inverse table with M_LPS_vs_M_LA

## The ebseq table is null.

## Working on 11/16: GM_LP_vs_GM_LA which is: GM_LP/GM_LA.

## Found table with GM_LP_vs_GM_LA

## The ebseq table is null.

## Working on 12/16: M_LP_vs_M_LA which is: M_LP/M_LA.

## Found table with M_LP_vs_M_LA

## The ebseq table is null.

## Working on 13/16: GM_NS_vs_M_NS which is: GM_NS/M_NS.

## Found inverse table with M_NS_vs_GM_NS

## The ebseq table is null.

## Working on 14/16: GM_LPS_vs_M_LPS which is: GM_LPS/M_LPS.

## Found inverse table with M_LPS_vs_GM_LPS

## The ebseq table is null.

## Working on 15/16: GM_LP_vs_M_LP which is: GM_LP/M_LP.

## Found inverse table with M_LP_vs_GM_LP

## The ebseq table is null.

## Working on 16/16: GM_LA_vs_M_LA which is: GM_LA/M_LA.

## Found inverse table with M_LA_vs_GM_LA

## The ebseq table is null.

## Adding venn plots for GM_LPS_vs_GM_NS.

## Limma expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LP_vs_GM_NS.

## Limma expression coefficients for GM_LP_vs_GM_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LP_vs_GM_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LP_vs_GM_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LA_vs_GM_NS.

## Limma expression coefficients for GM_LA_vs_GM_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LA_vs_GM_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LA_vs_GM_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LPS_vs_M_NS.

## Limma expression coefficients for M_LPS_vs_M_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LPS_vs_M_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LPS_vs_M_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LP_vs_M_NS.

## Limma expression coefficients for M_LP_vs_M_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LP_vs_M_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LP_vs_M_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LA_vs_M_NS.

## Limma expression coefficients for M_LA_vs_M_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LA_vs_M_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LA_vs_M_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LP_vs_GM_LPS.

## Limma expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LA_vs_GM_LPS.

## Limma expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LP_vs_M_LPS.

## Limma expression coefficients for M_LP_vs_M_LPS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LP_vs_M_LPS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LP_vs_M_LPS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LA_vs_M_LPS.

## Limma expression coefficients for M_LA_vs_M_LPS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LA_vs_M_LPS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LA_vs_M_LPS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LP_vs_GM_LA.

## Limma expression coefficients for GM_LP_vs_GM_LA; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LP_vs_GM_LA; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LP_vs_GM_LA; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for M_LP_vs_M_LA.

## Limma expression coefficients for M_LP_vs_M_LA; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for M_LP_vs_M_LA; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for M_LP_vs_M_LA; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_NS_vs_M_NS.

## Limma expression coefficients for GM_NS_vs_M_NS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_NS_vs_M_NS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_NS_vs_M_NS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LPS_vs_M_LPS.

## Limma expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LP_vs_M_LP.

## Limma expression coefficients for GM_LP_vs_M_LP; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LP_vs_M_LP; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LP_vs_M_LP; R^2: 0.996; equation: y = 0.997x + 0.0247

## Adding venn plots for GM_LA_vs_M_LA.

## Limma expression coefficients for GM_LA_vs_M_LA; R^2: 0.998; equation: y = 0.998x - 0.00439

## Deseq expression coefficients for GM_LA_vs_M_LA; R^2: 0.992; equation: y = 0.982x + 0.172

## Edger expression coefficients for GM_LA_vs_M_LA; R^2: 0.996; equation: y = 0.997x + 0.0247

## Writing summary information, compare_plot is: TRUE.

## Performing save of excel/pat_gs_sva_tables-v20200330.xlsx.

pat_gs_sva_sig <- extract_significant_genes(
    pat_gs_sva_tables,
    excel=glue::glue("excel/pat_gs_sva_sig-v{ver}.xlsx"))

## Writing a legend of columns.

## Did not find the ebseq_logfc, skipping ebseq.

## Writing excel data according to limma for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 313 genes.

## After (adj)p filter, the down genes table has 162 genes.

## After fold change filter, the up genes table has 190 genes.

## After fold change filter, the down genes table has 27 genes.

## Writing excel data according to limma for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 457 genes.

## After (adj)p filter, the down genes table has 379 genes.

## After fold change filter, the up genes table has 304 genes.

## After fold change filter, the down genes table has 80 genes.

## Writing excel data according to limma for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 522 genes.

## After (adj)p filter, the down genes table has 480 genes.

## After fold change filter, the up genes table has 283 genes.

## After fold change filter, the down genes table has 91 genes.

## Writing excel data according to limma for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 2033 genes.

## After (adj)p filter, the down genes table has 2942 genes.

## After fold change filter, the up genes table has 993 genes.

## After fold change filter, the down genes table has 560 genes.

## Writing excel data according to limma for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 2195 genes.

## After (adj)p filter, the down genes table has 2904 genes.

## After fold change filter, the up genes table has 1036 genes.

## After fold change filter, the down genes table has 712 genes.

## Writing excel data according to limma for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 2336 genes.

## After (adj)p filter, the down genes table has 2909 genes.

## After fold change filter, the up genes table has 1038 genes.

## After fold change filter, the down genes table has 765 genes.

## Writing excel data according to limma for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 122 genes.

## After (adj)p filter, the down genes table has 88 genes.

## After fold change filter, the up genes table has 73 genes.

## After fold change filter, the down genes table has 11 genes.

## Writing excel data according to limma for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 5 genes.

## After (adj)p filter, the down genes table has 2 genes.

## After fold change filter, the up genes table has 5 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to limma for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 1034 genes.

## After (adj)p filter, the down genes table has 1032 genes.

## After fold change filter, the up genes table has 253 genes.

## After fold change filter, the down genes table has 255 genes.

## Writing excel data according to limma for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 621 genes.

## After (adj)p filter, the down genes table has 513 genes.

## After fold change filter, the up genes table has 145 genes.

## After fold change filter, the down genes table has 132 genes.

## Writing excel data according to limma for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to limma for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 143 genes.

## After (adj)p filter, the down genes table has 162 genes.

## After fold change filter, the up genes table has 55 genes.

## After fold change filter, the down genes table has 28 genes.

## Writing excel data according to limma for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 1072 genes.

## After (adj)p filter, the down genes table has 1161 genes.

## After fold change filter, the up genes table has 482 genes.

## After fold change filter, the down genes table has 388 genes.

## Writing excel data according to limma for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 2105 genes.

## After (adj)p filter, the down genes table has 1894 genes.

## After fold change filter, the up genes table has 687 genes.

## After fold change filter, the down genes table has 757 genes.

## Writing excel data according to limma for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 1740 genes.

## After (adj)p filter, the down genes table has 1764 genes.

## After fold change filter, the up genes table has 712 genes.

## After fold change filter, the down genes table has 657 genes.

## Writing excel data according to limma for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 1933 genes.

## After (adj)p filter, the down genes table has 2036 genes.

## After fold change filter, the up genes table has 717 genes.

## After fold change filter, the down genes table has 701 genes.

## Printing significant genes to the file: excel/pat_gs_sva_sig-v20200330.xlsx

## 1/16: Creating significant table up_limma_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_limma_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_limma_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_limma_M_LPS_vs_M_NS

## 5/16: Creating significant table up_limma_M_LP_vs_M_NS

## 6/16: Creating significant table up_limma_M_LA_vs_M_NS

## 7/16: Creating significant table up_limma_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_limma_GM_LA_vs_GM_LPS

## The down table GM_LA_vs_GM_LPS is empty.

## 9/16: Creating significant table up_limma_M_LP_vs_M_LPS

## 10/16: Creating significant table up_limma_M_LA_vs_M_LPS

## The up table GM_LP_vs_GM_LA is empty.

## The down table GM_LP_vs_GM_LA is empty.

## 12/16: Creating significant table up_limma_M_LP_vs_M_LA

## 13/16: Creating significant table up_limma_GM_NS_vs_M_NS

## 14/16: Creating significant table up_limma_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_limma_GM_LP_vs_M_LP

## 16/16: Creating significant table up_limma_GM_LA_vs_M_LA

## Writing excel data according to edger for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 513 genes.

## After (adj)p filter, the down genes table has 155 genes.

## After fold change filter, the up genes table has 281 genes.

## After fold change filter, the down genes table has 34 genes.

## Writing excel data according to edger for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 672 genes.

## After (adj)p filter, the down genes table has 415 genes.

## After fold change filter, the up genes table has 420 genes.

## After fold change filter, the down genes table has 113 genes.

## Writing excel data according to edger for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 702 genes.

## After (adj)p filter, the down genes table has 447 genes.

## After fold change filter, the up genes table has 390 genes.

## After fold change filter, the down genes table has 127 genes.

## Writing excel data according to edger for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 2003 genes.

## After (adj)p filter, the down genes table has 1889 genes.

## After fold change filter, the up genes table has 1124 genes.

## After fold change filter, the down genes table has 583 genes.

## Writing excel data according to edger for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 2033 genes.

## After (adj)p filter, the down genes table has 2035 genes.

## After fold change filter, the up genes table has 1160 genes.

## After fold change filter, the down genes table has 675 genes.

## Writing excel data according to edger for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 2127 genes.

## After (adj)p filter, the down genes table has 2066 genes.

## After fold change filter, the up genes table has 1185 genes.

## After fold change filter, the down genes table has 741 genes.

## Writing excel data according to edger for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 219 genes.

## After (adj)p filter, the down genes table has 190 genes.

## After fold change filter, the up genes table has 136 genes.

## After fold change filter, the down genes table has 29 genes.

## Writing excel data according to edger for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 20 genes.

## After (adj)p filter, the down genes table has 12 genes.

## After fold change filter, the up genes table has 17 genes.

## After fold change filter, the down genes table has 7 genes.

## Writing excel data according to edger for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 837 genes.

## After (adj)p filter, the down genes table has 1134 genes.

## After fold change filter, the up genes table has 304 genes.

## After fold change filter, the down genes table has 299 genes.

## Writing excel data according to edger for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 513 genes.

## After (adj)p filter, the down genes table has 674 genes.

## After fold change filter, the up genes table has 201 genes.

## After fold change filter, the down genes table has 173 genes.

## Writing excel data according to edger for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 40 genes.

## After (adj)p filter, the down genes table has 8 genes.

## After fold change filter, the up genes table has 25 genes.

## After fold change filter, the down genes table has 2 genes.

## Writing excel data according to edger for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 136 genes.

## After (adj)p filter, the down genes table has 141 genes.

## After fold change filter, the up genes table has 52 genes.

## After fold change filter, the down genes table has 36 genes.

## Writing excel data according to edger for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 1192 genes.

## After (adj)p filter, the down genes table has 1013 genes.

## After fold change filter, the up genes table has 613 genes.

## After fold change filter, the down genes table has 419 genes.

## Writing excel data according to edger for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 1709 genes.

## After (adj)p filter, the down genes table has 1733 genes.

## After fold change filter, the up genes table has 741 genes.

## After fold change filter, the down genes table has 816 genes.

## Writing excel data according to edger for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 1543 genes.

## After (adj)p filter, the down genes table has 1487 genes.

## After fold change filter, the up genes table has 790 genes.

## After fold change filter, the down genes table has 698 genes.

## Writing excel data according to edger for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 1674 genes.

## After (adj)p filter, the down genes table has 1755 genes.

## After fold change filter, the up genes table has 793 genes.

## After fold change filter, the down genes table has 808 genes.

## Printing significant genes to the file: excel/pat_gs_sva_sig-v20200330.xlsx

## 1/16: Creating significant table up_edger_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_edger_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_edger_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_edger_M_LPS_vs_M_NS

## 5/16: Creating significant table up_edger_M_LP_vs_M_NS

## 6/16: Creating significant table up_edger_M_LA_vs_M_NS

## 7/16: Creating significant table up_edger_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_edger_GM_LA_vs_GM_LPS

## 9/16: Creating significant table up_edger_M_LP_vs_M_LPS

## 10/16: Creating significant table up_edger_M_LA_vs_M_LPS

## 11/16: Creating significant table up_edger_GM_LP_vs_GM_LA

## 12/16: Creating significant table up_edger_M_LP_vs_M_LA

## 13/16: Creating significant table up_edger_GM_NS_vs_M_NS

## 14/16: Creating significant table up_edger_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_edger_GM_LP_vs_M_LP

## 16/16: Creating significant table up_edger_GM_LA_vs_M_LA

## Writing excel data according to deseq for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 758 genes.

## After (adj)p filter, the down genes table has 333 genes.

## After fold change filter, the up genes table has 283 genes.

## After fold change filter, the down genes table has 26 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 986 genes.

## After (adj)p filter, the down genes table has 598 genes.

## After fold change filter, the up genes table has 437 genes.

## After fold change filter, the down genes table has 101 genes.

## Writing excel data according to deseq for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 1043 genes.

## After (adj)p filter, the down genes table has 693 genes.

## After fold change filter, the up genes table has 360 genes.

## After fold change filter, the down genes table has 106 genes.

## Writing excel data according to deseq for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 2627 genes.

## After (adj)p filter, the down genes table has 2685 genes.

## After fold change filter, the up genes table has 931 genes.

## After fold change filter, the down genes table has 455 genes.

## Writing excel data according to deseq for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 2664 genes.

## After (adj)p filter, the down genes table has 2712 genes.

## After fold change filter, the up genes table has 939 genes.

## After fold change filter, the down genes table has 530 genes.

## Writing excel data according to deseq for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 2782 genes.

## After (adj)p filter, the down genes table has 2811 genes.

## After fold change filter, the up genes table has 913 genes.

## After fold change filter, the down genes table has 587 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 350 genes.

## After (adj)p filter, the down genes table has 306 genes.

## After fold change filter, the up genes table has 140 genes.

## After fold change filter, the down genes table has 31 genes.

## Writing excel data according to deseq for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 24 genes.

## After (adj)p filter, the down genes table has 23 genes.

## After fold change filter, the up genes table has 14 genes.

## After fold change filter, the down genes table has 2 genes.

## Writing excel data according to deseq for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 1274 genes.

## After (adj)p filter, the down genes table has 1494 genes.

## After fold change filter, the up genes table has 254 genes.

## After fold change filter, the down genes table has 250 genes.

## Writing excel data according to deseq for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 778 genes.

## After (adj)p filter, the down genes table has 1006 genes.

## After fold change filter, the up genes table has 145 genes.

## After fold change filter, the down genes table has 155 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 55 genes.

## After (adj)p filter, the down genes table has 14 genes.

## After fold change filter, the up genes table has 34 genes.

## After fold change filter, the down genes table has 2 genes.

## Writing excel data according to deseq for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 309 genes.

## After (adj)p filter, the down genes table has 303 genes.

## After fold change filter, the up genes table has 57 genes.

## After fold change filter, the down genes table has 28 genes.

## Writing excel data according to deseq for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 1581 genes.

## After (adj)p filter, the down genes table has 1354 genes.

## After fold change filter, the up genes table has 519 genes.

## After fold change filter, the down genes table has 361 genes.

## Writing excel data according to deseq for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 2450 genes.

## After (adj)p filter, the down genes table has 2355 genes.

## After fold change filter, the up genes table has 615 genes.

## After fold change filter, the down genes table has 689 genes.

## Writing excel data according to deseq for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 2243 genes.

## After (adj)p filter, the down genes table has 2004 genes.

## After fold change filter, the up genes table has 665 genes.

## After fold change filter, the down genes table has 587 genes.

## Writing excel data according to deseq for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 2468 genes.

## After (adj)p filter, the down genes table has 2322 genes.

## After fold change filter, the up genes table has 677 genes.

## After fold change filter, the down genes table has 638 genes.

## Printing significant genes to the file: excel/pat_gs_sva_sig-v20200330.xlsx

## 1/16: Creating significant table up_deseq_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_deseq_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_deseq_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_deseq_M_LPS_vs_M_NS

## 5/16: Creating significant table up_deseq_M_LP_vs_M_NS

## 6/16: Creating significant table up_deseq_M_LA_vs_M_NS

## 7/16: Creating significant table up_deseq_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_deseq_GM_LA_vs_GM_LPS

## 9/16: Creating significant table up_deseq_M_LP_vs_M_LPS

## 10/16: Creating significant table up_deseq_M_LA_vs_M_LPS

## 11/16: Creating significant table up_deseq_GM_LP_vs_GM_LA

## 12/16: Creating significant table up_deseq_M_LP_vs_M_LA

## 13/16: Creating significant table up_deseq_GM_NS_vs_M_NS

## 14/16: Creating significant table up_deseq_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_deseq_GM_LP_vs_M_LP

## 16/16: Creating significant table up_deseq_GM_LA_vs_M_LA

## Writing excel data according to basic for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 1 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 1 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 1056 genes.

## After (adj)p filter, the down genes table has 1023 genes.

## After fold change filter, the up genes table has 623 genes.

## After fold change filter, the down genes table has 268 genes.

## Writing excel data according to basic for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 1090 genes.

## After (adj)p filter, the down genes table has 1024 genes.

## After fold change filter, the up genes table has 636 genes.

## After fold change filter, the down genes table has 357 genes.

## Writing excel data according to basic for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 1446 genes.

## After (adj)p filter, the down genes table has 1500 genes.

## After fold change filter, the up genes table has 695 genes.

## After fold change filter, the down genes table has 452 genes.

## Writing excel data according to basic for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 233 genes.

## After (adj)p filter, the down genes table has 222 genes.

## After fold change filter, the up genes table has 104 genes.

## After fold change filter, the down genes table has 85 genes.

## Writing excel data according to basic for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 2 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 1 genes.

## Writing excel data according to basic for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 1 genes.

## After (adj)p filter, the down genes table has 1 genes.

## After fold change filter, the up genes table has 1 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 36 genes.

## After (adj)p filter, the down genes table has 34 genes.

## After fold change filter, the up genes table has 26 genes.

## After fold change filter, the down genes table has 16 genes.

## Writing excel data according to basic for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 163 genes.

## After (adj)p filter, the down genes table has 262 genes.

## After fold change filter, the up genes table has 98 genes.

## After fold change filter, the down genes table has 188 genes.

## Writing excel data according to basic for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 230 genes.

## After (adj)p filter, the down genes table has 311 genes.

## After fold change filter, the up genes table has 141 genes.

## After fold change filter, the down genes table has 198 genes.

## Writing excel data according to basic for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 289 genes.

## After (adj)p filter, the down genes table has 378 genes.

## After fold change filter, the up genes table has 174 genes.

## After fold change filter, the down genes table has 235 genes.

## Printing significant genes to the file: excel/pat_gs_sva_sig-v20200330.xlsx

## 1/16: Creating significant table up_basic_GM_LPS_vs_GM_NS

## The down table GM_LPS_vs_GM_NS is empty.

## The up table GM_LP_vs_GM_NS is empty.

## The down table GM_LP_vs_GM_NS is empty.

## The up table GM_LA_vs_GM_NS is empty.

## The down table GM_LA_vs_GM_NS is empty.

## 4/16: Creating significant table up_basic_M_LPS_vs_M_NS

## 5/16: Creating significant table up_basic_M_LP_vs_M_NS

## 6/16: Creating significant table up_basic_M_LA_vs_M_NS

## The up table GM_LP_vs_GM_LPS is empty.

## The down table GM_LP_vs_GM_LPS is empty.

## The up table GM_LA_vs_GM_LPS is empty.

## The down table GM_LA_vs_GM_LPS is empty.

## 9/16: Creating significant table up_basic_M_LP_vs_M_LPS

## The up table M_LA_vs_M_LPS is empty.

## The up table GM_LP_vs_GM_LA is empty.

## The down table GM_LP_vs_GM_LA is empty.

## 12/16: Creating significant table up_basic_M_LP_vs_M_LA

## The down table M_LP_vs_M_LA is empty.

## 13/16: Creating significant table up_basic_GM_NS_vs_M_NS

## 14/16: Creating significant table up_basic_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_basic_GM_LP_vs_M_LP

## 16/16: Creating significant table up_basic_GM_LA_vs_M_LA

## Adding significance bar plots.

6.0.2 Patient as batch, growth and stimulation as condition, batch in model.

pat_gs_batch <- set_expt_conditions(hs_filt, fact="gp")
pat_gs_batch <- set_expt_batches(pat_gs_batch, fact="patient")
pat_gs_batch_de <- all_pairwise(pat_gs_batch, model_batch=TRUE)

## Plotting a PCA before surrogates/batch inclusion.

## Not putting labels on the plot.

## Using limma's removeBatchEffect to visualize with(out) batch inclusion.

## Not putting labels on the plot.

## Finished running DE analyses, collecting outputs.

## Comparing analyses.

pat_gs_batch_tables <- combine_de_tables(
    pat_gs_batch_de, keepers=keepers,
    excel=glue::glue("excel/pat_gs_batch_tables-v{ver}.xlsx"))

## Deleting the file excel/pat_gs_batch_tables-v20200330.xlsx before writing the tables.

## Writing a legend of columns.

## Printing a pca plot before/after surrogates/batch estimation.

## Working on 1/16: GM_LPS_vs_GM_NS which is: GM_LPS/GM_NS.

## Found inverse table with GM_NS_vs_GM_LPS

## The ebseq table is null.

## Working on 2/16: GM_LP_vs_GM_NS which is: GM_LP/GM_NS.

## Found inverse table with GM_NS_vs_GM_LP

## The ebseq table is null.

## Working on 3/16: GM_LA_vs_GM_NS which is: GM_LA/GM_NS.

## Found inverse table with GM_NS_vs_GM_LA

## The ebseq table is null.

## Working on 4/16: M_LPS_vs_M_NS which is: M_LPS/M_NS.

## Found inverse table with M_NS_vs_M_LPS

## The ebseq table is null.

## Working on 5/16: M_LP_vs_M_NS which is: M_LP/M_NS.

## Found inverse table with M_NS_vs_M_LP

## The ebseq table is null.

## Working on 6/16: M_LA_vs_M_NS which is: M_LA/M_NS.

## Found inverse table with M_NS_vs_M_LA

## The ebseq table is null.

## Working on 7/16: GM_LP_vs_GM_LPS which is: GM_LP/GM_LPS.

## Found inverse table with GM_LPS_vs_GM_LP

## The ebseq table is null.

## Working on 8/16: GM_LA_vs_GM_LPS which is: GM_LA/GM_LPS.

## Found inverse table with GM_LPS_vs_GM_LA

## The ebseq table is null.

## Working on 9/16: M_LP_vs_M_LPS which is: M_LP/M_LPS.

## Found inverse table with M_LPS_vs_M_LP

## The ebseq table is null.

## Working on 10/16: M_LA_vs_M_LPS which is: M_LA/M_LPS.

## Found inverse table with M_LPS_vs_M_LA

## The ebseq table is null.

## Working on 11/16: GM_LP_vs_GM_LA which is: GM_LP/GM_LA.

## Found table with GM_LP_vs_GM_LA

## The ebseq table is null.

## Working on 12/16: M_LP_vs_M_LA which is: M_LP/M_LA.

## Found table with M_LP_vs_M_LA

## The ebseq table is null.

## Working on 13/16: GM_NS_vs_M_NS which is: GM_NS/M_NS.

## Found inverse table with M_NS_vs_GM_NS

## The ebseq table is null.

## Working on 14/16: GM_LPS_vs_M_LPS which is: GM_LPS/M_LPS.

## Found inverse table with M_LPS_vs_GM_LPS

## The ebseq table is null.

## Working on 15/16: GM_LP_vs_M_LP which is: GM_LP/M_LP.

## Found inverse table with M_LP_vs_GM_LP

## The ebseq table is null.

## Working on 16/16: GM_LA_vs_M_LA which is: GM_LA/M_LA.

## Found inverse table with M_LA_vs_GM_LA

## The ebseq table is null.

## Adding venn plots for GM_LPS_vs_GM_NS.

## Limma expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LPS_vs_GM_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LP_vs_GM_NS.

## Limma expression coefficients for GM_LP_vs_GM_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LP_vs_GM_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LP_vs_GM_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LA_vs_GM_NS.

## Limma expression coefficients for GM_LA_vs_GM_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LA_vs_GM_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LA_vs_GM_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LPS_vs_M_NS.

## Limma expression coefficients for M_LPS_vs_M_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LPS_vs_M_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LPS_vs_M_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LP_vs_M_NS.

## Limma expression coefficients for M_LP_vs_M_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LP_vs_M_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LP_vs_M_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LA_vs_M_NS.

## Limma expression coefficients for M_LA_vs_M_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LA_vs_M_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LA_vs_M_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LP_vs_GM_LPS.

## Limma expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LP_vs_GM_LPS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LA_vs_GM_LPS.

## Limma expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LA_vs_GM_LPS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LP_vs_M_LPS.

## Limma expression coefficients for M_LP_vs_M_LPS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LP_vs_M_LPS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LP_vs_M_LPS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LA_vs_M_LPS.

## Limma expression coefficients for M_LA_vs_M_LPS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LA_vs_M_LPS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LA_vs_M_LPS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LP_vs_GM_LA.

## Limma expression coefficients for GM_LP_vs_GM_LA; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LP_vs_GM_LA; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LP_vs_GM_LA; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for M_LP_vs_M_LA.

## Limma expression coefficients for M_LP_vs_M_LA; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for M_LP_vs_M_LA; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for M_LP_vs_M_LA; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_NS_vs_M_NS.

## Limma expression coefficients for GM_NS_vs_M_NS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_NS_vs_M_NS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_NS_vs_M_NS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LPS_vs_M_LPS.

## Limma expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LPS_vs_M_LPS; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LP_vs_M_LP.

## Limma expression coefficients for GM_LP_vs_M_LP; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LP_vs_M_LP; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LP_vs_M_LP; R^2: 0.997; equation: y = 1x - 0.0147

## Adding venn plots for GM_LA_vs_M_LA.

## Limma expression coefficients for GM_LA_vs_M_LA; R^2: 0.998; equation: y = 0.999x - 0.00812

## Deseq expression coefficients for GM_LA_vs_M_LA; R^2: 0.99; equation: y = 1x - 0.0224

## Edger expression coefficients for GM_LA_vs_M_LA; R^2: 0.997; equation: y = 1x - 0.0147

## Writing summary information, compare_plot is: TRUE.

## Performing save of excel/pat_gs_batch_tables-v20200330.xlsx.

pat_gs_batch_sig <- extract_significant_genes(
    pat_gs_batch_tables,
    excel=glue::glue("excel/pat_gs_batch_sig-v{ver}.xlsx"))

## Writing a legend of columns.

## Did not find the ebseq_logfc, skipping ebseq.

## Writing excel data according to limma for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 300 genes.

## After (adj)p filter, the down genes table has 123 genes.

## After fold change filter, the up genes table has 199 genes.

## After fold change filter, the down genes table has 27 genes.

## Writing excel data according to limma for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 565 genes.

## After (adj)p filter, the down genes table has 538 genes.

## After fold change filter, the up genes table has 365 genes.

## After fold change filter, the down genes table has 143 genes.

## Writing excel data according to limma for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 472 genes.

## After (adj)p filter, the down genes table has 313 genes.

## After fold change filter, the up genes table has 285 genes.

## After fold change filter, the down genes table has 80 genes.

## Writing excel data according to limma for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 1867 genes.

## After (adj)p filter, the down genes table has 2465 genes.

## After fold change filter, the up genes table has 992 genes.

## After fold change filter, the down genes table has 563 genes.

## Writing excel data according to limma for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 1935 genes.

## After (adj)p filter, the down genes table has 2362 genes.

## After fold change filter, the up genes table has 1017 genes.

## After fold change filter, the down genes table has 682 genes.

## Writing excel data according to limma for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 2073 genes.

## After (adj)p filter, the down genes table has 2511 genes.

## After fold change filter, the up genes table has 1030 genes.

## After fold change filter, the down genes table has 778 genes.

## Writing excel data according to limma for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 116 genes.

## After (adj)p filter, the down genes table has 71 genes.

## After fold change filter, the up genes table has 76 genes.

## After fold change filter, the down genes table has 14 genes.

## Writing excel data according to limma for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 3 genes.

## After (adj)p filter, the down genes table has 1 genes.

## After fold change filter, the up genes table has 3 genes.

## After fold change filter, the down genes table has 1 genes.

## Writing excel data according to limma for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 817 genes.

## After (adj)p filter, the down genes table has 837 genes.

## After fold change filter, the up genes table has 261 genes.

## After fold change filter, the down genes table has 249 genes.

## Writing excel data according to limma for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 385 genes.

## After (adj)p filter, the down genes table has 335 genes.

## After fold change filter, the up genes table has 129 genes.

## After fold change filter, the down genes table has 117 genes.

## Writing excel data according to limma for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 2 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 2 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to limma for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 95 genes.

## After (adj)p filter, the down genes table has 91 genes.

## After fold change filter, the up genes table has 38 genes.

## After fold change filter, the down genes table has 25 genes.

## Writing excel data according to limma for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 864 genes.

## After (adj)p filter, the down genes table has 901 genes.

## After fold change filter, the up genes table has 440 genes.

## After fold change filter, the down genes table has 365 genes.

## Writing excel data according to limma for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 1817 genes.

## After (adj)p filter, the down genes table has 1670 genes.

## After fold change filter, the up genes table has 709 genes.

## After fold change filter, the down genes table has 779 genes.

## Writing excel data according to limma for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 2050 genes.

## After (adj)p filter, the down genes table has 1963 genes.

## After fold change filter, the up genes table has 847 genes.

## After fold change filter, the down genes table has 818 genes.

## Writing excel data according to limma for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 1887 genes.

## After (adj)p filter, the down genes table has 1817 genes.

## After fold change filter, the up genes table has 809 genes.

## After fold change filter, the down genes table has 742 genes.

## Printing significant genes to the file: excel/pat_gs_batch_sig-v20200330.xlsx

## 1/16: Creating significant table up_limma_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_limma_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_limma_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_limma_M_LPS_vs_M_NS

## 5/16: Creating significant table up_limma_M_LP_vs_M_NS

## 6/16: Creating significant table up_limma_M_LA_vs_M_NS

## 7/16: Creating significant table up_limma_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_limma_GM_LA_vs_GM_LPS

## 9/16: Creating significant table up_limma_M_LP_vs_M_LPS

## 10/16: Creating significant table up_limma_M_LA_vs_M_LPS

## 11/16: Creating significant table up_limma_GM_LP_vs_GM_LA

## The down table GM_LP_vs_GM_LA is empty.

## 12/16: Creating significant table up_limma_M_LP_vs_M_LA

## 13/16: Creating significant table up_limma_GM_NS_vs_M_NS

## 14/16: Creating significant table up_limma_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_limma_GM_LP_vs_M_LP

## 16/16: Creating significant table up_limma_GM_LA_vs_M_LA

## Writing excel data according to edger for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 414 genes.

## After (adj)p filter, the down genes table has 102 genes.

## After fold change filter, the up genes table has 267 genes.

## After fold change filter, the down genes table has 31 genes.

## Writing excel data according to edger for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 700 genes.

## After (adj)p filter, the down genes table has 484 genes.

## After fold change filter, the up genes table has 431 genes.

## After fold change filter, the down genes table has 168 genes.

## Writing excel data according to edger for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 600 genes.

## After (adj)p filter, the down genes table has 291 genes.

## After fold change filter, the up genes table has 377 genes.

## After fold change filter, the down genes table has 106 genes.

## Writing excel data according to edger for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 1969 genes.

## After (adj)p filter, the down genes table has 1647 genes.

## After fold change filter, the up genes table has 1153 genes.

## After fold change filter, the down genes table has 533 genes.

## Writing excel data according to edger for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 1909 genes.

## After (adj)p filter, the down genes table has 1688 genes.

## After fold change filter, the up genes table has 1164 genes.

## After fold change filter, the down genes table has 608 genes.

## Writing excel data according to edger for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 1973 genes.

## After (adj)p filter, the down genes table has 1896 genes.

## After fold change filter, the up genes table has 1170 genes.

## After fold change filter, the down genes table has 723 genes.

## Writing excel data according to edger for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 189 genes.

## After (adj)p filter, the down genes table has 133 genes.

## After fold change filter, the up genes table has 126 genes.

## After fold change filter, the down genes table has 37 genes.

## Writing excel data according to edger for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 18 genes.

## After (adj)p filter, the down genes table has 4 genes.

## After fold change filter, the up genes table has 15 genes.

## After fold change filter, the down genes table has 2 genes.

## Writing excel data according to edger for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 667 genes.

## After (adj)p filter, the down genes table has 855 genes.

## After fold change filter, the up genes table has 283 genes.

## After fold change filter, the down genes table has 284 genes.

## Writing excel data according to edger for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 313 genes.

## After (adj)p filter, the down genes table has 443 genes.

## After fold change filter, the up genes table has 164 genes.

## After fold change filter, the down genes table has 151 genes.

## Writing excel data according to edger for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 48 genes.

## After (adj)p filter, the down genes table has 14 genes.

## After fold change filter, the up genes table has 36 genes.

## After fold change filter, the down genes table has 7 genes.

## Writing excel data according to edger for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 91 genes.

## After (adj)p filter, the down genes table has 69 genes.

## After fold change filter, the up genes table has 45 genes.

## After fold change filter, the down genes table has 23 genes.

## Writing excel data according to edger for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 900 genes.

## After (adj)p filter, the down genes table has 831 genes.

## After fold change filter, the up genes table has 558 genes.

## After fold change filter, the down genes table has 406 genes.

## Writing excel data according to edger for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 1470 genes.

## After (adj)p filter, the down genes table has 1702 genes.

## After fold change filter, the up genes table has 720 genes.

## After fold change filter, the down genes table has 890 genes.

## Writing excel data according to edger for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 1722 genes.

## After (adj)p filter, the down genes table has 1849 genes.

## After fold change filter, the up genes table has 865 genes.

## After fold change filter, the down genes table has 935 genes.

## Writing excel data according to edger for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 1682 genes.

## After (adj)p filter, the down genes table has 1612 genes.

## After fold change filter, the up genes table has 835 genes.

## After fold change filter, the down genes table has 827 genes.

## Printing significant genes to the file: excel/pat_gs_batch_sig-v20200330.xlsx

## 1/16: Creating significant table up_edger_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_edger_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_edger_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_edger_M_LPS_vs_M_NS

## 5/16: Creating significant table up_edger_M_LP_vs_M_NS

## 6/16: Creating significant table up_edger_M_LA_vs_M_NS

## 7/16: Creating significant table up_edger_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_edger_GM_LA_vs_GM_LPS

## 9/16: Creating significant table up_edger_M_LP_vs_M_LPS

## 10/16: Creating significant table up_edger_M_LA_vs_M_LPS

## 11/16: Creating significant table up_edger_GM_LP_vs_GM_LA

## 12/16: Creating significant table up_edger_M_LP_vs_M_LA

## 13/16: Creating significant table up_edger_GM_NS_vs_M_NS

## 14/16: Creating significant table up_edger_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_edger_GM_LP_vs_M_LP

## 16/16: Creating significant table up_edger_GM_LA_vs_M_LA

## Writing excel data according to deseq for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 625 genes.

## After (adj)p filter, the down genes table has 269 genes.

## After fold change filter, the up genes table has 283 genes.

## After fold change filter, the down genes table has 31 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 1305 genes.

## After (adj)p filter, the down genes table has 870 genes.

## After fold change filter, the up genes table has 442 genes.

## After fold change filter, the down genes table has 182 genes.

## Writing excel data according to deseq for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 906 genes.

## After (adj)p filter, the down genes table has 543 genes.

## After fold change filter, the up genes table has 367 genes.

## After fold change filter, the down genes table has 107 genes.

## Writing excel data according to deseq for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 2396 genes.

## After (adj)p filter, the down genes table has 2147 genes.

## After fold change filter, the up genes table has 976 genes.

## After fold change filter, the down genes table has 419 genes.

## Writing excel data according to deseq for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 2308 genes.

## After (adj)p filter, the down genes table has 2184 genes.

## After fold change filter, the up genes table has 950 genes.

## After fold change filter, the down genes table has 507 genes.

## Writing excel data according to deseq for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 2396 genes.

## After (adj)p filter, the down genes table has 2402 genes.

## After fold change filter, the up genes table has 931 genes.

## After fold change filter, the down genes table has 589 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 275 genes.

## After (adj)p filter, the down genes table has 228 genes.

## After fold change filter, the up genes table has 125 genes.

## After fold change filter, the down genes table has 31 genes.

## Writing excel data according to deseq for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 16 genes.

## After (adj)p filter, the down genes table has 14 genes.

## After fold change filter, the up genes table has 10 genes.

## After fold change filter, the down genes table has 4 genes.

## Writing excel data according to deseq for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 984 genes.

## After (adj)p filter, the down genes table has 1280 genes.

## After fold change filter, the up genes table has 237 genes.

## After fold change filter, the down genes table has 269 genes.

## Writing excel data according to deseq for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 491 genes.

## After (adj)p filter, the down genes table has 738 genes.

## After fold change filter, the up genes table has 131 genes.

## After fold change filter, the down genes table has 153 genes.

## Writing excel data according to deseq for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 72 genes.

## After (adj)p filter, the down genes table has 17 genes.

## After fold change filter, the up genes table has 41 genes.

## After fold change filter, the down genes table has 2 genes.

## Writing excel data according to deseq for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 191 genes.

## After (adj)p filter, the down genes table has 161 genes.

## After fold change filter, the up genes table has 51 genes.

## After fold change filter, the down genes table has 24 genes.

## Writing excel data according to deseq for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 1208 genes.

## After (adj)p filter, the down genes table has 1073 genes.

## After fold change filter, the up genes table has 509 genes.

## After fold change filter, the down genes table has 362 genes.

## Writing excel data according to deseq for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 2075 genes.

## After (adj)p filter, the down genes table has 2112 genes.

## After fold change filter, the up genes table has 644 genes.

## After fold change filter, the down genes table has 749 genes.

## Writing excel data according to deseq for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 2477 genes.

## After (adj)p filter, the down genes table has 2236 genes.

## After fold change filter, the up genes table has 742 genes.

## After fold change filter, the down genes table has 731 genes.

## Writing excel data according to deseq for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 2232 genes.

## After (adj)p filter, the down genes table has 2014 genes.

## After fold change filter, the up genes table has 725 genes.

## After fold change filter, the down genes table has 687 genes.

## Printing significant genes to the file: excel/pat_gs_batch_sig-v20200330.xlsx

## 1/16: Creating significant table up_deseq_GM_LPS_vs_GM_NS

## 2/16: Creating significant table up_deseq_GM_LP_vs_GM_NS

## 3/16: Creating significant table up_deseq_GM_LA_vs_GM_NS

## 4/16: Creating significant table up_deseq_M_LPS_vs_M_NS

## 5/16: Creating significant table up_deseq_M_LP_vs_M_NS

## 6/16: Creating significant table up_deseq_M_LA_vs_M_NS

## 7/16: Creating significant table up_deseq_GM_LP_vs_GM_LPS

## 8/16: Creating significant table up_deseq_GM_LA_vs_GM_LPS

## 9/16: Creating significant table up_deseq_M_LP_vs_M_LPS

## 10/16: Creating significant table up_deseq_M_LA_vs_M_LPS

## 11/16: Creating significant table up_deseq_GM_LP_vs_GM_LA

## 12/16: Creating significant table up_deseq_M_LP_vs_M_LA

## 13/16: Creating significant table up_deseq_GM_NS_vs_M_NS

## 14/16: Creating significant table up_deseq_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_deseq_GM_LP_vs_M_LP

## 16/16: Creating significant table up_deseq_GM_LA_vs_M_LA

## Writing excel data according to basic for GM_LPS_vs_GM_NS: 1/64.

## After (adj)p filter, the up genes table has 1 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 1 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LP_vs_GM_NS: 2/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LA_vs_GM_NS: 3/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LPS_vs_M_NS: 4/64.

## After (adj)p filter, the up genes table has 1056 genes.

## After (adj)p filter, the down genes table has 1023 genes.

## After fold change filter, the up genes table has 623 genes.

## After fold change filter, the down genes table has 268 genes.

## Writing excel data according to basic for M_LP_vs_M_NS: 5/64.

## After (adj)p filter, the up genes table has 1090 genes.

## After (adj)p filter, the down genes table has 1024 genes.

## After fold change filter, the up genes table has 636 genes.

## After fold change filter, the down genes table has 357 genes.

## Writing excel data according to basic for M_LA_vs_M_NS: 6/64.

## After (adj)p filter, the up genes table has 1446 genes.

## After (adj)p filter, the down genes table has 1500 genes.

## After fold change filter, the up genes table has 695 genes.

## After fold change filter, the down genes table has 452 genes.

## Writing excel data according to basic for GM_LP_vs_GM_LPS: 7/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_LA_vs_GM_LPS: 8/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LP_vs_M_LPS: 9/64.

## After (adj)p filter, the up genes table has 233 genes.

## After (adj)p filter, the down genes table has 222 genes.

## After fold change filter, the up genes table has 104 genes.

## After fold change filter, the down genes table has 85 genes.

## Writing excel data according to basic for M_LA_vs_M_LPS: 10/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 2 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 1 genes.

## Writing excel data according to basic for GM_LP_vs_GM_LA: 11/64.

## After (adj)p filter, the up genes table has 0 genes.

## After (adj)p filter, the down genes table has 0 genes.

## After fold change filter, the up genes table has 0 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for M_LP_vs_M_LA: 12/64.

## After (adj)p filter, the up genes table has 1 genes.

## After (adj)p filter, the down genes table has 1 genes.

## After fold change filter, the up genes table has 1 genes.

## After fold change filter, the down genes table has 0 genes.

## Writing excel data according to basic for GM_NS_vs_M_NS: 13/64.

## After (adj)p filter, the up genes table has 36 genes.

## After (adj)p filter, the down genes table has 34 genes.

## After fold change filter, the up genes table has 26 genes.

## After fold change filter, the down genes table has 16 genes.

## Writing excel data according to basic for GM_LPS_vs_M_LPS: 14/64.

## After (adj)p filter, the up genes table has 163 genes.

## After (adj)p filter, the down genes table has 262 genes.

## After fold change filter, the up genes table has 98 genes.

## After fold change filter, the down genes table has 188 genes.

## Writing excel data according to basic for GM_LP_vs_M_LP: 15/64.

## After (adj)p filter, the up genes table has 230 genes.

## After (adj)p filter, the down genes table has 311 genes.

## After fold change filter, the up genes table has 141 genes.

## After fold change filter, the down genes table has 198 genes.

## Writing excel data according to basic for GM_LA_vs_M_LA: 16/64.

## After (adj)p filter, the up genes table has 289 genes.

## After (adj)p filter, the down genes table has 378 genes.

## After fold change filter, the up genes table has 174 genes.

## After fold change filter, the down genes table has 235 genes.

## Printing significant genes to the file: excel/pat_gs_batch_sig-v20200330.xlsx

## 1/16: Creating significant table up_basic_GM_LPS_vs_GM_NS

## The down table GM_LPS_vs_GM_NS is empty.

## The up table GM_LP_vs_GM_NS is empty.

## The down table GM_LP_vs_GM_NS is empty.

## The up table GM_LA_vs_GM_NS is empty.

## The down table GM_LA_vs_GM_NS is empty.

## 4/16: Creating significant table up_basic_M_LPS_vs_M_NS

## 5/16: Creating significant table up_basic_M_LP_vs_M_NS

## 6/16: Creating significant table up_basic_M_LA_vs_M_NS

## The up table GM_LP_vs_GM_LPS is empty.

## The down table GM_LP_vs_GM_LPS is empty.

## The up table GM_LA_vs_GM_LPS is empty.

## The down table GM_LA_vs_GM_LPS is empty.

## 9/16: Creating significant table up_basic_M_LP_vs_M_LPS

## The up table M_LA_vs_M_LPS is empty.

## The up table GM_LP_vs_GM_LA is empty.

## The down table GM_LP_vs_GM_LA is empty.

## 12/16: Creating significant table up_basic_M_LP_vs_M_LA

## The down table M_LP_vs_M_LA is empty.

## 13/16: Creating significant table up_basic_GM_NS_vs_M_NS

## 14/16: Creating significant table up_basic_GM_LPS_vs_M_LPS

## 15/16: Creating significant table up_basic_GM_LP_vs_M_LP

## 16/16: Creating significant table up_basic_GM_LA_vs_M_LA

## Adding significance bar plots.

6.1 Compare sva/batch in model

comp <- compare_de_results(pat_gs_sva_tables, pat_gs_batch_tables)

## Testing method: limma.

## Adding method: limma to the set.

## Testing method: deseq.

## Adding method: deseq to the set.

## Testing method: edger.

## Adding method: edger to the set.

## Testing method: ebseq.

## The first datum is missing method: ebseq.

## Testing method: basic.

## Adding method: basic to the set.

##  Starting method limma, table GM_LPS_vs_GM_NS.

##  Starting method limma, table GM_LP_vs_GM_NS.

##  Starting method limma, table GM_LA_vs_GM_NS.

##  Starting method limma, table M_LPS_vs_M_NS.

##  Starting method limma, table M_LP_vs_M_NS.

##  Starting method limma, table M_LA_vs_M_NS.

##  Starting method limma, table GM_LP_vs_GM_LPS.

##  Starting method limma, table GM_LA_vs_GM_LPS.

##  Starting method limma, table M_LP_vs_M_LPS.

##  Starting method limma, table M_LA_vs_M_LPS.

##  Starting method limma, table GM_LP_vs_GM_LA.

##  Starting method limma, table M_LP_vs_M_LA.

##  Starting method limma, table GM_NS_vs_M_NS.

##  Starting method limma, table GM_LPS_vs_M_LPS.

##  Starting method limma, table GM_LP_vs_M_LP.

##  Starting method limma, table GM_LA_vs_M_LA.

##  Starting method deseq, table GM_LPS_vs_GM_NS.

##  Starting method deseq, table GM_LP_vs_GM_NS.

##  Starting method deseq, table GM_LA_vs_GM_NS.

##  Starting method deseq, table M_LPS_vs_M_NS.

##  Starting method deseq, table M_LP_vs_M_NS.

##  Starting method deseq, table M_LA_vs_M_NS.

##  Starting method deseq, table GM_LP_vs_GM_LPS.

##  Starting method deseq, table GM_LA_vs_GM_LPS.

##  Starting method deseq, table M_LP_vs_M_LPS.

##  Starting method deseq, table M_LA_vs_M_LPS.

##  Starting method deseq, table GM_LP_vs_GM_LA.

##  Starting method deseq, table M_LP_vs_M_LA.

##  Starting method deseq, table GM_NS_vs_M_NS.

##  Starting method deseq, table GM_LPS_vs_M_LPS.

##  Starting method deseq, table GM_LP_vs_M_LP.

##  Starting method deseq, table GM_LA_vs_M_LA.

##  Starting method edger, table GM_LPS_vs_GM_NS.

##  Starting method edger, table GM_LP_vs_GM_NS.

##  Starting method edger, table GM_LA_vs_GM_NS.

##  Starting method edger, table M_LPS_vs_M_NS.

##  Starting method edger, table M_LP_vs_M_NS.

##  Starting method edger, table M_LA_vs_M_NS.

##  Starting method edger, table GM_LP_vs_GM_LPS.

##  Starting method edger, table GM_LA_vs_GM_LPS.

##  Starting method edger, table M_LP_vs_M_LPS.

##  Starting method edger, table M_LA_vs_M_LPS.

##  Starting method edger, table GM_LP_vs_GM_LA.

##  Starting method edger, table M_LP_vs_M_LA.

##  Starting method edger, table GM_NS_vs_M_NS.

##  Starting method edger, table GM_LPS_vs_M_LPS.

##  Starting method edger, table GM_LP_vs_M_LP.

##  Starting method edger, table GM_LA_vs_M_LA.

##  Starting method basic, table GM_LPS_vs_GM_NS.

##  Starting method basic, table GM_LP_vs_GM_NS.

##  Starting method basic, table GM_LA_vs_GM_NS.

##  Starting method basic, table M_LPS_vs_M_NS.

##  Starting method basic, table M_LP_vs_M_NS.

##  Starting method basic, table M_LA_vs_M_NS.

##  Starting method basic, table GM_LP_vs_GM_LPS.

##  Starting method basic, table GM_LA_vs_GM_LPS.

##  Starting method basic, table M_LP_vs_M_LPS.

##  Starting method basic, table M_LA_vs_M_LPS.

##  Starting method basic, table GM_LP_vs_GM_LA.

##  Starting method basic, table M_LP_vs_M_LA.

##  Starting method basic, table GM_NS_vs_M_NS.

##  Starting method basic, table GM_LPS_vs_M_LPS.

##  Starting method basic, table GM_LP_vs_M_LP.

##  Starting method basic, table GM_LA_vs_M_LA.

## Look at the logFC comparisons:
comp$lfc_heat

## Look at the p-value comparisons:
comp$p_heat

## It appears edgeR is a bit more sensitive to changes in the model.

7 Perform some ontology searches

It appears that I crashed the gProfiler web server by sending in my various searches. So I will leave these off for the moment and replace them with some clusterProfiler searches.

7.1 GM Comparisons

7.1.1 GM LPS vs GM NS

table <- "GM_LPS_vs_GM_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sa_sig[["deseq"]][["downs"]][[table]]

## Error in eval(expr, envir, enclos): object 'pat_gs_sa_sig' not found

de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 272 out of 283 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 26 MF, 1076 BP, and 2 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 565 enriched hits.

## Found 268 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 59 KEGG enriched hits.

## Attempting DAVID search.

## Loading required package: RDAVIDWebService

## Loading required package: graph

## Loading required package: GOstats

## Loading required package: Category

## Loading required package: stats4

## Loading required package: AnnotationDbi

## Loading required package: IRanges

## Loading required package: S4Vectors

## 
## Attaching package: 'S4Vectors'

## The following object is masked from 'package:base':
## 
##     expand.grid

## Loading required package: Matrix

## 
## Attaching package: 'Matrix'

## The following object is masked from 'package:S4Vectors':
## 
##     expand

## 
## Attaching package: 'GOstats'

## The following object is masked from 'package:AnnotationDbi':
## 
##     makeGOGraph

## Loading required package: ggplot2

## Need help? Try Stackoverflow: https://stackoverflow.com/tags/ggplot2

## 
## Attaching package: 'RDAVIDWebService'

## The following object is masked from 'package:AnnotationDbi':
## 
##     species

## The following object is masked from 'package:IRanges':
## 
##     members

## The following objects are masked from 'package:BiocGenerics':
## 
##     counts, species, type

## Found 550 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Error in rownames(sig_genes): object 'down' not found

## No significant hits going down
ont_down$plots$ego_all_mf

## Error in eval(expr, envir, enclos): object 'ont_down' not found

ont_down$plots$ego_all_bp

## Error in eval(expr, envir, enclos): object 'ont_down' not found

ont_down$plots$dot_all_mf

## Error in eval(expr, envir, enclos): object 'ont_down' not found

ont_down$plots$dot_all_bp

## Error in eval(expr, envir, enclos): object 'ont_down' not found

ont_down$plots$map_all_mf

## Error in eval(expr, envir, enclos): object 'ont_down' not found

ont_down$plots$map_all_bp

## Error in eval(expr, envir, enclos): object 'ont_down' not found

7.1.2 GM LP vs GM NS

table <- "GM_LP_vs_GM_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 423 out of 437 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 22 MF, 1070 BP, and 1 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 524 enriched hits.

## Found 414 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 63 KEGG enriched hits.

## Attempting DAVID search.

## Found 566 DAVID hits.

## Plotting results.

## Using `stress` as default layout

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 98 out of 101 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 0 MF, 0 BP, and 0 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 746 enriched hits.

## Found 98 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 0 DAVID hits.

## Plotting results.

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.1.3 GM LA vs GM NS

table <- "GM_LA_vs_GM_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 341 out of 360 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 20 MF, 977 BP, and 1 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 489 enriched hits.

## Found 334 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 63 KEGG enriched hits.

## Attempting DAVID search.

## Found 503 DAVID hits.

## Plotting results.

## Using `stress` as default layout

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 103 out of 106 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 2 MF, 0 BP, and 0 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 824 enriched hits.

## Found 102 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 0 DAVID hits.

## Plotting results.

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.2 GM Comparisons against LPS

7.2.1 GM LP vs GM LPS

table <- "GM_LP_vs_GM_LPS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 135 out of 140 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 9 MF, 73 BP, and 1 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 133 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 1 KEGG enriched hits.

## Attempting DAVID search.

## Found 2 DAVID hits.

## Plotting results.

## Using `stress` as default layout

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 29 out of 31 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 12 MF, 76 BP, and 2 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 28 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 1 KEGG enriched hits.

## Attempting DAVID search.

## Found 0 DAVID hits.

## Plotting results.

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.2.2 GM LA vs GM LPS

table <- "GM_LA_vs_GM_LPS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
## No hits!

7.2.3 GM LA vs GM LP

up <- pat_gs_sva_sig[["deseq"]][["ups"]][["GM_LP_vs_GM_LA"]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][["GM_LP_vs_GM_LA"]]

ont_up <- simple_gprofiler(up)

## Performing gProfiler GO search of 34 genes against hsapiens.

## Performing gProfiler KEGG search of 34 genes against hsapiens.

## Performing gProfiler REAC search of 34 genes against hsapiens.

## Performing gProfiler MI search of 34 genes against hsapiens.

## Performing gProfiler TF search of 34 genes against hsapiens.

## Performing gProfiler CORUM search of 34 genes against hsapiens.

## Performing gProfiler HP search of 34 genes against hsapiens.

ont_up[["pvalue_plots"]][["mfp_plot_over"]]

## Error in ont_up[["pvalue_plots"]][["mfp_plot_over"]]: subscript out of bounds

ont_up[["pvalue_plots"]][["bpp_plot_over"]]

## Error in ont_up[["pvalue_plots"]][["bpp_plot_over"]]: subscript out of bounds

ont_up[["pvalue_plots"]][["kegg_plot_over"]]

## Error in ont_up[["pvalue_plots"]][["kegg_plot_over"]]: subscript out of bounds

ont_up[["pvalue_plots"]][["reactome_plot_over"]]

## Error in ont_up[["pvalue_plots"]][["reactome_plot_over"]]: subscript out of bounds

ont_up[["pvalue_plots"]][["hp_plot_over"]]

## Error in ont_up[["pvalue_plots"]][["hp_plot_over"]]: subscript out of bounds

## Down only had 2 genes.

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 31 out of 34 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 20 MF, 41 BP, and 0 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 31 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 1 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 1 out of 2 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 6 MF, 144 BP, and 1 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 1 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 32 KEGG enriched hits.

## Attempting DAVID search.

## Warning in clusterProfiler::enrichDAVID(gene = sig_gene_list, minGSSize =
## min_groupsize, : No significant enrichment found...

## Found 0 DAVID hits.

## Plotting results.

## Using `stress` as default layout

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.3 M Comparisons

7.3.1 M LPS vs M NS

This search crashed the gProfiler server, so I will stop it at least for the moment.

table <- "M_LPS_vs_M_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 908 out of 931 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 41 MF, 1455 BP, and 12 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 834 enriched hits.

## Found 891 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 83 KEGG enriched hits.

## Attempting DAVID search.

## Found 824 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 449 out of 455 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 4 MF, 5 BP, and 0 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 861 enriched hits.

## Found 441 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 0 DAVID hits.

## Plotting results.

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.3.2 M LP vs M NS

table <- "M_LP_vs_M_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 913 out of 939 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 38 MF, 1267 BP, and 5 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 507 enriched hits.

## Found 898 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 74 KEGG enriched hits.

## Attempting DAVID search.

## Found 720 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 524 out of 530 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 16 MF, 88 BP, and 8 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 517 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 10 DAVID hits.

## Plotting results.

## No significant hits going down
ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.3.3 M LA vs M NS

table <- "M_LA_vs_M_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_gprofiler(up)
ont_up[["pvalue_plots"]][["mfp_plot_over"]]
ont_up[["pvalue_plots"]][["bpp_plot_over"]]
ont_up[["pvalue_plots"]][["kegg_plot_over"]]
ont_up[["pvalue_plots"]][["reactome_plot_over"]]
ont_up[["pvalue_plots"]][["hp_plot_over"]]

ont_down <- simple_gprofiler(down)
ont_down[["pvalue_plots"]][["mfp_plot_over"]]
ont_down[["pvalue_plots"]][["bpp_plot_over"]]
ont_down[["pvalue_plots"]][["kegg_plot_over"]]
ont_down[["pvalue_plots"]][["reactome_plot_over"]]
ont_down[["pvalue_plots"]][["hp_plot_over"]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 892 out of 913 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 25 MF, 1394 BP, and 4 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 592 enriched hits.

## Found 873 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 80 KEGG enriched hits.

## Attempting DAVID search.

## Found 763 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 584 out of 587 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 18 MF, 51 BP, and 10 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 360 enriched hits.

## Found 575 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 15 DAVID hits.

## Plotting results.

ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

7.4 M vs GM

table <- "GM_NS_vs_M_NS"
up <- pat_gs_sva_sig[["deseq"]][["ups"]][[table]]
down <- pat_gs_sva_sig[["deseq"]][["downs"]][[table]]
de <- pat_gs_sva_tables[["data"]][[table]]

ont_up <- simple_clusterprofiler(sig_genes=up, de_table=de,
                                 do_david=TRUE, david_user="abelew@umd.edu",
                                 fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 501 out of 519 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 26 MF, 374 BP, and 10 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 488 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 8 KEGG enriched hits.

## Attempting DAVID search.

## Found 133 DAVID hits.

## Plotting results.

ont_up$plots$ego_sig_mf

ont_up$plots$ego_sig_bp

ont_up$plots$dot_sig_mf

ont_up$plots$dot_sig_bp

ont_up$plots$map_sig_mf

ont_up$plots$map_sig_bp

ont_down <- simple_clusterprofiler(sig_genes=down, de_table=de,
                                   do_david=TRUE, david_user="abelew@umd.edu",
                                   fc_column="deseq_logfc", orgdb="org.Hs.eg.db")

## Testing available OrgDb keytypes for the best mapping to entrez.

## Chose keytype: ENSEMBL for all genes because it had 20051 out of 34431 genes.

## Chose keytype: ENSEMBL for sig genes because it had 338 out of 361 genes.

## Calculating GO groups.

## Found 146 MF, 558 BP, and 755 CC hits.

## Calculating enriched GO groups.

## Found 9 MF, 92 BP, and 2 CC enriched hits.

## Performing GSE analyses of gene lists (this is slow).

## Found 0 enriched hits.

## Found 332 matches between the significant gene list and kegg universe.

## Performing KEGG analyses.

## Found 19329 matches between the gene list and kegg universe.

## Found 0 KEGG enriched hits.

## Attempting DAVID search.

## Found 47 DAVID hits.

## Plotting results.

ont_down$plots$ego_all_mf

ont_down$plots$ego_all_bp

ont_down$plots$dot_all_mf

ont_down$plots$dot_all_bp

ont_down$plots$map_all_mf

## NULL

ont_down$plots$map_all_bp

## NULL

pander::pander(sessionInfo())
message(paste0("This is hpgltools commit: ", get_git_commit()))
this_save <- paste0(gsub(pattern="\\.Rmd", replace="", x=rmd_file), "-v", ver, ".rda.xz")
message(paste0("Saving to ", this_save))
tmp <- sm(saveme(filename=this_save))

loadme(filename=this_save)