Tuesday, September 1, 2015

9 1 2015 Mean Ct Value Stats and Graphs R script

Steven has been helping me cull and curate the Ct values from the qPCR runs over the summer. He's developed a spreadsheet of mean Ct values of the good replicate runs. I've taken those values and created delta Ct values of the target over Actin. Then I log transformed the data to produce as normal data as possible. Finally I ran ANOVA and Tukey's Honest Significant Difference Post Hoc, on the log transformed data. Then created boxplots of the delta Ct values. I will create a table of the significant differences as a clear way to present the results. Below is the annotated code for the R script which I'll post a link to in the github repository at the bottom.

#Necessary Packages to manipulate data and plot values.
require(plyr)
require(ggplot2)
require(splitstackshape)

#Read in mean Ct value table
dCt<-read.csv("CTvalues83115.csv", header=T)
#Split SAMPLE_ID column to create columns for population, treatment, and sample number
dCt<-cSplit(dCt,"SAMPLE_ID", sep= "_", drop=F)
#rename columns appropriately
dCt<-rename(dCt,replace=c("SAMPLE_ID_1"="Pop","SAMPLE_ID_2"="Treat","SAMPLE_ID_3"="Sample"))

#divide each target of interest by the mean Ct value of the Actin Normalizing gene
dCt$CARM<-(dCt$CarmmeanCt/dCt$Actinmeanct)
dCt$TLR<-(dCt$TLRaverage/dCt$Actinmeanct)
dCt$CRAF<-(dCt$CRAFctaverage/dCt$Actinmeanct)
dCt$H2AV<-(dCt$H2AVavgct/dCt$Actinmeanct)
dCt$PGRP<-(dCt$PGRPaverage/dCt$Actinmeanct)
dCt$HSP70<-(dCt$HSP70averageCt/dCt$Actinmeanct)
dCt$BMP2<-(dCt$BMP2average/dCt$Actinmeanct)
dCt$GRB2<-(dCt$GRB2average/dCt$Actinmeanct)
dCt$PGEEP4<-(dCt$PGEEP4ctav/dCt$Actinmeanct)

#log transform the data to develop normality in data
dCt$CARMlog<-log(dCt$CARM)
dCt$TLRlog<-log(dCt$TLR)
dCt$H2AVlog<-log(dCt$H2AV)
dCt$PGRPlog<-log(dCt$PGRP)
dCt$HSP70log<-log(dCt$HSP70)
dCt$BMP2log<-log(dCt$BMP2)
dCt$GRB2log<-log(dCt$GRB2)
dCt$PGEEP4log<-log(dCt$PGEEP4)
dCt$CRAFlog<-log(dCt$CRAF)

#Run ANOVA's on all log transformed data as well as Tukey's Honestly Significant Difference post hoc test
CARM<-aov(CARMlog~Pop+Treat+Pop:Treat, data=dCt)
CARM
TukeyHSD(CARM)

TLR<-aov(TLRlog~Pop+Treat+Pop:Treat, data=dCt)
TLR
TukeyHSD(TLR)

H2AV<-aov(H2AVlog~Pop+Treat+Pop:Treat, data=dCt)
H2AV
TukeyHSD(H2AV)

PGRP<-aov(PGRPlog~Pop+Treat+Pop:Treat, data=dCt)
PGRP
TukeyHSD(PGRP)

HSP70<-aov(HSP70log~Pop+Treat+Pop:Treat, data=dCt)
HSP70
TukeyHSD(HSP70)

BMP2<-aov(BMP2log~Pop+Treat+Pop:Treat, data=dCt)
BMP2
TukeyHSD(BMP2)

GRB2<-aov(GRB2log~Pop+Treat+Pop:Treat, data=dCt)
GRB2
TukeyHSD(GRB2)

PGEEP4<-aov(PGEEP4log~Pop+Treat+Pop:Treat, data=dCt)
PGEEP4
TukeyHSD(PGEEP4)

CRAF<-aov(CRAFlog~Pop+Treat+Pop:Treat, data=dCt)
CRAF
TukeyHSD(CRAF)

#graph all raw mean Ct values to produce boxplots to visualize data

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=CARM,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")


ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=TLR, fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=H2AV,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=PGRP,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=HSP70,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=BMP2,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=GRB2,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=PGEEP4,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

ggplot(data=dCt)+geom_boxplot(aes(x=Treat, y=CRAF,fill=Pop))+theme_bw()+
  theme(axis.text.x=element_text(size=20), axis.text.y=element_text(size=20),
        axis.title.x=element_text(size=25), axis.title.y=element_text(size=25),
        legend.position=c(.1,.1),panel.grid.major=element_blank())+
  labs(x="Treatment", y="target/actin delta Ct")

Graphs produced










Significant Differences Table
CARMlogTLRlogH2AVlogPGRPlogHSP70logBMP2logGRB2logPGEEP4logCRAFlog
Overall
C:TXXX
C:MXXX
T:MXXXXX
N:HX
H:SX
S:NX
Control
N:H
H:S
S:N
Temperature
N:H
H:S
S:N
Mechanical
N:H
H:S
S:NX
Control:Temp
N
H
SX
Control:Mech
N
H
SXX
Temp:Mech
N
H
SXXX

You can find the raw mean Ct values here and the R script here.

September 2015 Goals

Its September 1st, Time for GOALLLLLLLLSSSSSSSSS!!!!! Let's take a quick look at the goals for last month and what I accomplished of them.

August 2015 Goals

1.Complete Analysis of qPCR data
2. Write Results
3. Complete Chapter 2
4. Complete Thesis
5. Schedule New Defense Date
6. Get TA position for the fall
7. Drink some delicious Flying Bike beer. 

How I did:

1. I'm really close. We've got most if not all of the qPCR data that we're going to use for chapter 2. I have written an R script to run through the average Ct data and produce statistics and graphs for it. I'm tidying it up today so I can basically say this one is done. 

2. This will be done in the next two weeks now that the Analysis is basically done. 

3. This will be done in the next two weeks.

4. See 2 & 3.

5. New Defense date is October 14th at 11 am.

6. Still trying to find one. Waiting to hear about the Intro R TAship.

7. The brewery opened and the beer is delicious!


September 2015 Goals

1.Complete Results, Chapter 2, Thesis
2. Sending Draft of thesis to committee (mid Sept)
3. Begin preparations for defense presentation
4. Find TAship/Start TAing
5. Attend PCSGA (Network, view other presentations, work on presentation)
6. Submit Chapter 1 to JSR
7. Go camping!