I am trying to make a gene co-expression network using RNA-seq data to study the host-pathogen interaction. My paired end RNA seq data consist of three time points - Day1(2 replicates), Day3(2 replicates) and Day 7(2 replicates)

I processed the RNA-seq reads and aligned it to the corresponding genome. I used HTseq-count to get the read counts and was able to pick the differentially expressed genes using DEseq2 R package. My DEseq2 output files (Day1 vs Day3) and (Day1 vs Day7) where Day1 is taken as control and Day 3 and Day 7 is taken as treated. This analysis returned 13,780 gene as DEGs and 16,530 genes as DEGs after filtering using q-value cutoff < 0.1.

I am trying construct the gene expression matrix and co-expression network using the function cor() in R language.

My queries are as follows

1. How to construct the gene expression matrix.
2. What should be taken as input to construct Gene Expression matrix

I am following the methodology from the paper published in Frontiers of Genetics, November 2016, volume 7, article 206, " A Network Approach of Gene Co-expression in the Zea mays/Aspergillus flavus Pathosystem to Map Host/Pathogen Interaction Pathways".