Search
Question: (Closed) FW: EdgeR GLM
0
6.5 years ago by
Shona Wood70
Shona Wood70 wrote:
Hi, I am trying to look at gene expression changes with increasing age. I have three time points young, middle age and old. I have been using the following, ordered factors and coefficient in order to get a fold change and FDR for a linear change in gene expression with increasing age. Can you look over my code and see if I have done this right? Is age.L the fold change in gene expression with increasing age? > library(edgeR) > library(limma) > targets<-read.delim (file="coding_targets.txt") > targets$age<-factor(ordered (targets$age)) > targets$sample<-factor(targets$sample) > targets X files age sample 1 p1 p1_coding_counts.txt a-young 1 2 p2 p2_coding_counts.txt a-young 2 3 p3 p3_coding_counts.txt a-young 3 4 p7 p7_coding_counts.txt b-middle 7 5 p8 p8_coding_counts.txt b-middle 8 6 p9 p9_coding_counts.txt b-middle 9 7 p4 p4_coding_counts.txt c-old 4 8 p5 p5_coding_counts.txt c-old 5 9 p6 p6_coding_counts.txt c-old 6 > d<-readDGE(targets, skip=1, comment.char='#') > colnames(d)<-row.names(targets) > d<- calcNormFactors(d) > d An object of class "DGEList" $samples X files age sample lib.size norm.factors 1 p1 p1_coding_counts.txt a-young 1 3445622 0.9655724 2 p2 p2_coding_counts.txt a-young 2 2696547 0.9902573 3 p3 p3_coding_counts.txt a-young 3 3308099 0.9787044 4 p7 p7_coding_counts.txt b-middle 7 2503479 1.0584660 5 p8 p8_coding_counts.txt b-middle 8 2639127 1.0477893 6 p9 p9_coding_counts.txt b-middle 9 2696547 0.9902573 7 p4 p4_coding_counts.txt c-old 4 3037440 0.9778553 8 p5 p5_coding_counts.txt c-old 5 2647915 1.0144109 9 p6 p6_coding_counts.txt c-old 6 2475370 0.9809077$counts 1 2 3 4 5 6 7 8 9 ENSRNOG00000000007 1287 1285 1041 788 968 1285 1092 1009 960 ENSRNOG00000000008 0 0 0 0 1 0 0 1 0 ENSRNOG00000000009 0 0 0 3 0 0 1 0 0 ENSRNOG00000000010 38 405 50 18 105 405 372 42 282 ENSRNOG00000000012 0 0 0 0 0 0 0 0 0 22932 more rows ... > d<-d[rowSums(d$counts)>9,] > design<- model.matrix(~ age, data = targets) > design (Intercept) age.L age.Q 1 1 -7.071068e-01 0.4082483 2 1 -7.071068e-01 0.4082483 3 1 -7.071068e-01 0.4082483 4 1 -7.850462e-17 -0.8164966 5 1 -7.850462e-17 -0.8164966 6 1 -7.850462e-17 -0.8164966 7 1 7.071068e-01 0.4082483 8 1 7.071068e-01 0.4082483 9 1 7.071068e-01 0.4082483 attr(,"assign") [1] 0 1 1 attr(,"contrasts") attr(,"contrasts")$age [1] "contr.poly" > d<-estimateGLMCommonDisp(d, design) > names(d) [1] "samples" "counts" "common.dispersion" > glmfit<- glmFit(d, design, dispersion=d\$common.dispersion) > results<- glmLRT(d, glmfit, coef=c(2,3)) > topTags(results) Coefficient: age.L age.Q logConc age.L age.Q LR P.Value ENSRNOG00000011672 -11.27122 -3.2445886 -1.9642177 109.53641 1.638591e-24 Many Thanks [[alternative HTML version deleted]]
written 6.5 years ago by Shona Wood70

Hello Shona Wood!

Questions similar to yours can already be found at:

We have closed your question to allow us to keep similar content in the same thread.

If you disagree with this please tell us why in a reply below. We'll be happy to talk about it.

Cheers!