Currently just using the following setup to look at differences in treatment across all 3 time points I have (30, 180, 360 min) however would like to also break down analysis by time, any suggestions for design would be much appreciated thank you!

hyMData<-read.csv("hyM.csv")

rownames(hyMData) <- hyMData[,1]

hyMData[,1] <- NULL

sampleNames<-c("30_V1","30_V2","30_V3","30_E1","30_E2","30_E3","180_V1","180_V2","180_V3","180_E1","180_E2","180_E3","360_V1","360_V2","360_V3","360_E1","360_E2","360_E3")
sampleTime<-c("30","30","30","30","30","30","180","180","180","180","180","180","360","360","360","360","360","360")
sampleBatch<-c("1","2","3","1","2","3","1","2","3","1","2","3","1","2","3","1","2","3")
sampleTreatment<-c("CO","CO","CO","EV","EV","EV","CO","CO","CO","EV","EV","EV","CO","CO","CO","EV","EV","EV")
colData<-data.frame(sampleName=sampleNames, time=sampleTime, hormone=sampleTreatment, batch=sampleBatch)

dds<-DESeqDataSetFromMatrix(countData = hyMData,
                            colData = colData,
                            design = ~ batch + treatment)
dds$treatment<-relevel(dds$treatment, ref="CO")
dds<-DESeq(dds, betaPrior=TRUE)
treatment<-results(dds, contrast=c("treatment","EV","CO"))
write.csv(treatment,"hyM.csv")

I think you might as well start by making a PCA plot of your samples...if you find that batch isn't having an effect, you might drop it as a design factor.

Take a look at the time series example in the workflow, it sounds similar to what you have:

https://bioconductor.org/packages/rnaseqGene