Question: limma complex contrast matrix
0
gravatar for Prasad Siddavatam
8.4 years ago by
United States
Prasad Siddavatam150 wrote:
Dear Mailing list, I Have couple of questions about contrast matrix I have 206 arrays and here are the sample targets looks like num sample age gender strain stage dis(time) --- ------- --- ------ ------ ------ ----- 1 GSM614020 12 Male DENV1 ucd fu 2 GSM614021 7 Male DENV1 ucd acute 4 GSM614023 12 Male DENV1 dss fu 5 GSM614024 11 Male DENV1 ucd dis 6 GSM614025 5 Male DENV1 dss 01 7 GSM614026 5 Male DENV1 dss fu 8 GSM614027 11 Male DENV2 ucd acute 9 GSM614028 4 Male DENV1 ucd dis 10 GSM614029 10 Female DENV1 dss acute 11 GSM614030 10 Female DENV1 dss 01 12 GSM614031 10 Female DENV1 dss dis 13 GSM614032 10 Female DENV1 dss fu 15 GSM614034 14 Male DENV1 dss fu 16 GSM614035 10 Female DENV1 dss dis 17 GSM614036 10 Female DENV1 dss 01 18 GSM614037 29 Female DENV1 dss acute 19 GSM614038 29 Female DENV1 dss 01 20 GSM614039 29 Female DENV1 dss dis 21 GSM614040 14 Male DENV1 dss acute 61 GSM614080 11 Male DENV2 dss acute 62 GSM614081 11 Male DENV2 dss dis 63 GSM614082 19 Male DENV1 dss acute This gives me 2 strain factors, 2 stage factors and 4 dis(time) factors. Among four time factors every time point is to be compared with "fu" time because this sample is taken in the follow up visit after recovery. (16 comparisons) I want to have DE genes from many contrasts 1. Expressed in DENV1 only and DEN2 only 2. Expressed in DENV1 with ucd and DEN1 with dss (Same for DEN2) 3. Expressed in DENV1 with ucd and acute (and for other times also) 4. Expressed in and so on...... TS <- paste(target$strain, target$stage, target$dis, sep = "."); TS <- factor(TS, levels = unique(TS)); design <- model.matrix(~0+TS); colnames(design) <- levels(TS); Here are the unique(TS) values are......... DENV1.ucd.fu DENV1.ucd.acute DENV1.ucd.dis DENV1.dss.fu DENV1.dss.01 DENV2.ucd.acute DENV1.dss.acute DENV1.dss.dis DENV2.ucd.fu DENV2.dss.acute DENV2.dss.fu DENV1.ucd.01 DENV2.ucd.01 DENV2.dss.dis DENV2.ucd.dis Here is my confusing "makeContrasts" function contrast <- makeContrasts(d1.dss.dis = "DENV1.dss.dis - DENV1.dss.fu", d1.dss.acute = "DENV1.dss.acute - DENV1.dss.fu", d1.dss.01 = "DENV1.dss.01 - DENV1.dss.fu", d1.ucd.dis = "DENV1.ucd.dis - DENV1.ucd.fu", d1.ucd.acute = "DENV1.ucd.acute - DENV1.ucd.fu", d1.ucd.01 = "DENV1.ucd.01 - DENV1.ucd.fu", d1.dssVSucd = "((DENV1.dss.acute + DENV1.dss.dis + DENV1.dss.01)/3 - (DENV1.dss.fu)) - ((DENV1.ucd.acute + DENV1.ucd.dis + DENV1.ucd.01)/3 - (DENV1.ucd.fu))", d2.dss.dis = "DENV2.dss.dis - DENV2.dss.fu", d2.dss.acute = "DENV2.dss.acute - DENV2.dss.fu", d2.ucd.dis = "DENV2.ucd.dis - DENV2.ucd.fu", d2.ucd.acute = "DENV2.ucd.acute - DENV2.ucd.fu", d2.ucd.01 = "DENV2.ucd.01 - DENV2.ucd.fu", d2.dssVSucd = "((DENV2.dss.acute + DENV2.dss.dis)/2 - (DENV2.dss.fu)) - ((DENV2.ucd.acute + DENV2.ucd.dis + DENV2.ucd.01)/3 - (DENV2.ucd.fu))", d1VSd2 = "(((DENV1.dss.acute + DENV1.dss.dis + DENV1.dss.01)/3 - (DENV1.dss.fu)) + ((DENV1.ucd.acute + DENV1.ucd.dis + DENV1.ucd.01)/3 - (DENV1.ucd.fu))) - ((DENV2.dss.acute + DENV2.dss.dis)/2 - (DENV2.dss.fu)) + ((DENV2.ucd.acute + DENV2.ucd.dis + DENV2.ucd.01)/3 - (DENV2.ucd.fu)))", levels = design); here is the actual contrasts matrix I got (columns 4 to 7 are pasted) Levels d1.ucd.dis d1.ucd.acute d1.ucd.01 d1.dssVSucd DENV1.ucd.fu -1 -1 -1 1.0000000 DENV1.ucd.acute 0 1 0 -0.3333333 DENV1.ucd.dis 1 0 0 -0.3333333 DENV1.dss.fu 0 0 0 -1.0000000 DENV1.dss.01 0 0 0 0.3333333 DENV2.ucd.acute 0 0 0 0.0000000 DENV1.dss.acute 0 0 0 0.3333333 DENV1.dss.dis 0 0 0 0.3333333 DENV2.ucd.fu 0 0 0 0.0000000 DENV2.dss.acute 0 0 0 0.0000000 DENV2.dss.fu 0 0 0 0.0000000 DENV1.ucd.01 0 0 1 -0.3333333 DENV2.ucd.01 0 0 0 0.0000000 DENV2.dss.dis 0 0 0 0.0000000 DENV2.ucd.dis 0 0 0 0.0000000 Please advise me whether I am doing it correct or not? Any help is highly appreciated Above all, I want to add "age" and "gender" factors also later. I can't even imagine how messy the contrast would become. Technically (statistically) can I add those? Thank you very much Prasad
dss • 600 views
ADD COMMENTlink modified 8.4 years ago by Ina Hoeschele610 • written 8.4 years ago by Prasad Siddavatam150
Answer: limma complex contrast matrix
0
gravatar for Ina Hoeschele
8.4 years ago by
Ina Hoeschele610
United States
Ina Hoeschele610 wrote:
Hi Gordon et al., I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? Many thanks, Ina
ADD COMMENTlink written 8.4 years ago by Ina Hoeschele610
Dear Ina: You do not need to specify positive control probes. They were automatically read in when you use read.ilmn function to read in both regular probe data and control probe data. You can use the code below to see the numbers of regular probes, negative control probes and positive control probes: x <- read.ilmn(files="probe profile.txt",ctrlfiles="control probe profile.txt", other.columns="Detection") # you will need to changes file names to the names of your files table(x$genes$Status) The neqc function will then perform a quantile normalization for data included in x using all the probes (including regular probes, negative controls and positive controls). The two parameters negctrl and regular tell neqc what the identifiers of regular probes and negative control probes are. The negative control probes play a particular role in the background correction (normexp background correction aided by negative controls). By default, regular probes have an identifier of "regular". Negative control probes by default have an identifier of "negative" which was the identifier used by almost all the BeadChip data we have ever seen. The neqc function should work in most cases with its default setting. For more details, please refer to the limma user guide for analyzing Illumina BeadChip data (section 11.7). Hope this helps. Cheers, Wei On Mar 23, 2011, at 6:40 AM, Ina Hoeschele wrote: > Hi Gordon et al., > I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? > Many thanks, Ina > > _______________________________________________ > Bioconductor mailing list > Bioconductor at r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:6}}
ADD REPLYlink written 8.4 years ago by Wei Shi3.1k
Hi Wei, thanks much for your quick response. However, I should have mentioned that I am not reading the data using read.ilmn, I am reading and pre-processing the bead-level data via beadarray (doing QC etc.), then summarizing the data and then calling neqc. Ina ----- Original Message ----- From: "Wei Shi" <shi@wehi.edu.au> To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> Cc: bioconductor at stat.math.ethz.ch Sent: Tuesday, March 22, 2011 6:08:34 PM Subject: Re: [BioC] limma neqc Dear Ina: You do not need to specify positive control probes. They were automatically read in when you use read.ilmn function to read in both regular probe data and control probe data. You can use the code below to see the numbers of regular probes, negative control probes and positive control probes: x <- read.ilmn(files="probe profile.txt",ctrlfiles="control probe profile.txt", other.columns="Detection") # you will need to changes file names to the names of your files table(x$genes$Status) The neqc function will then perform a quantile normalization for data included in x using all the probes (including regular probes, negative controls and positive controls). The two parameters negctrl and regular tell neqc what the identifiers of regular probes and negative control probes are. The negative control probes play a particular role in the background correction (normexp background correction aided by negative controls). By default, regular probes have an identifier of "regular". Negative control probes by default have an identifier of "negative" which was the identifier used by almost all the BeadChip data we have ever seen. The neqc function should work in most cases with its default setting. For more details, please refer to the limma user guide for analyzing Illumina BeadChip data (section 11.7). Hope this helps. Cheers, Wei On Mar 23, 2011, at 6:40 AM, Ina Hoeschele wrote: > Hi Gordon et al., > I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? > Many thanks, Ina > > _______________________________________________ > Bioconductor mailing list > Bioconductor at r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
ADD REPLYlink written 8.4 years ago by Ina Hoeschele610
Hi Ina: You will have to provide the "status" parameter to neqc then. This is a character vector which tells neqc which probes are regular one, which are negative control and which are positive controls. Regular probes should have a value of "regular" (same with value of the "regular" parameter) and negative controls should have a value of "negative". For positive controls, you can use their original names like "Housekeeping", "Biotin" etc, or you can just use a value of "positive" for all of them. The output of neqc will be a matrix which only includes regular probes. Please let me know if you have any further questions. Cheers, Wei On Mar 24, 2011, at 1:15 AM, Ina Hoeschele wrote: > Hi Wei, > thanks much for your quick response. However, I should have mentioned that I am not reading the data using read.ilmn, I am reading and pre-processing the bead-level data via beadarray (doing QC etc.), then summarizing the data and then calling neqc. > Ina > > ----- Original Message ----- > From: "Wei Shi" <shi at="" wehi.edu.au=""> > To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> > Cc: bioconductor at stat.math.ethz.ch > Sent: Tuesday, March 22, 2011 6:08:34 PM > Subject: Re: [BioC] limma neqc > > Dear Ina: > > You do not need to specify positive control probes. They were automatically read in when you use read.ilmn function to read in both regular probe data and control probe data. You can use the code below to see the numbers of regular probes, negative control probes and positive control probes: > > x <- read.ilmn(files="probe profile.txt",ctrlfiles="control probe profile.txt", other.columns="Detection") # you will need to changes file names to the names of your files > table(x$genes$Status) > > The neqc function will then perform a quantile normalization for data included in x using all the probes (including regular probes, negative controls and positive controls). The two parameters negctrl and regular tell neqc what the identifiers of regular probes and negative control probes are. The negative control probes play a particular role in the background correction (normexp background correction aided by negative controls). By default, regular probes have an identifier of "regular". Negative control probes by default have an identifier of "negative" which was the identifier used by almost all the BeadChip data we have ever seen. The neqc function should work in most cases with its default setting. > > For more details, please refer to the limma user guide for analyzing Illumina BeadChip data (section 11.7). > > Hope this helps. > > Cheers, > Wei > > On Mar 23, 2011, at 6:40 AM, Ina Hoeschele wrote: > >> Hi Gordon et al., >> I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? >> Many thanks, Ina >> >> _______________________________________________ >> Bioconductor mailing list >> Bioconductor at r-project.org >> https://stat.ethz.ch/mailman/listinfo/bioconductor >> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:14}}
ADD REPLYlink written 8.4 years ago by Wei Shi3.1k
Hi Wei, thanks again for your response. I am not sure I understand it correctly: I do provide a status vector and I tell neqc that negative control probes are "negative". But I am not sure what regular means. There are positive control probes (biotin, cy3, ERCC,...) and gene or non-control probes which are denoted by "Gene". So do I have to set regular = "Gene"? Ina ----- Original Message ----- From: "Wei Shi" <shi@wehi.edu.au> To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> Cc: bioconductor at stat.math.ethz.ch Sent: Wednesday, March 23, 2011 6:14:50 PM Subject: Re: [BioC] limma neqc Hi Ina: You will have to provide the "status" parameter to neqc then. This is a character vector which tells neqc which probes are regular one, which are negative control and which are positive controls. Regular probes should have a value of "regular" (same with value of the "regular" parameter) and negative controls should have a value of "negative". For positive controls, you can use their original names like "Housekeeping", "Biotin" etc, or you can just use a value of "positive" for all of them. The output of neqc will be a matrix which only includes regular probes. Please let me know if you have any further questions. Cheers, Wei On Mar 24, 2011, at 1:15 AM, Ina Hoeschele wrote: > Hi Wei, > thanks much for your quick response. However, I should have mentioned that I am not reading the data using read.ilmn, I am reading and pre-processing the bead-level data via beadarray (doing QC etc.), then summarizing the data and then calling neqc. > Ina > > ----- Original Message ----- > From: "Wei Shi" <shi at="" wehi.edu.au=""> > To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> > Cc: bioconductor at stat.math.ethz.ch > Sent: Tuesday, March 22, 2011 6:08:34 PM > Subject: Re: [BioC] limma neqc > > Dear Ina: > > You do not need to specify positive control probes. They were automatically read in when you use read.ilmn function to read in both regular probe data and control probe data. You can use the code below to see the numbers of regular probes, negative control probes and positive control probes: > > x <- read.ilmn(files="probe profile.txt",ctrlfiles="control probe profile.txt", other.columns="Detection") # you will need to changes file names to the names of your files > table(x$genes$Status) > > The neqc function will then perform a quantile normalization for data included in x using all the probes (including regular probes, negative controls and positive controls). The two parameters negctrl and regular tell neqc what the identifiers of regular probes and negative control probes are. The negative control probes play a particular role in the background correction (normexp background correction aided by negative controls). By default, regular probes have an identifier of "regular". Negative control probes by default have an identifier of "negative" which was the identifier used by almost all the BeadChip data we have ever seen. The neqc function should work in most cases with its default setting. > > For more details, please refer to the limma user guide for analyzing Illumina BeadChip data (section 11.7). > > Hope this helps. > > Cheers, > Wei > > On Mar 23, 2011, at 6:40 AM, Ina Hoeschele wrote: > >> Hi Gordon et al., >> I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? >> Many thanks, Ina >> >> _______________________________________________ >> Bioconductor mailing list >> Bioconductor at r-project.org >> https://stat.ethz.ch/mailman/listinfo/bioconductor >> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:10}}
ADD REPLYlink written 8.4 years ago by Ina Hoeschele610
Hi Ina: "regular" probes are those functional probes on the array (~48,000 such probes in WG6 or HT12), which are neither negative controls nor positive controls. The reason why we use "regular" as the default status name of functional probes in neqc function is because the read.ilmn function (BeadChip data input function in limma) sets the status of functional probes to "regular". The status names of negative control probes and positive control probes given by read.ilmn function are the same as their names in the control profile file. However because you provide a matrix to neqc function, it does not really matter what the status names are as long as the three neqc parameters including "status", "negctrl" and "regular" are consistent. For example, if you use "Gene" as the status name of regular probes in your "status" vector you should set "regular" parameter to "Gene" as well. Hope this clarifies this issue. Cheers, Wei On Mar 31, 2011, at 6:58 AM, Ina Hoeschele wrote: > Hi Wei, > thanks again for your response. I am not sure I understand it correctly: I do provide a status vector and I tell neqc that negative control probes are "negative". But I am not sure what regular means. There are positive control probes (biotin, cy3, ERCC,...) and gene or non-control probes which are denoted by "Gene". So do I have to set regular = "Gene"? > Ina > > > ----- Original Message ----- > From: "Wei Shi" <shi at="" wehi.edu.au=""> > To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> > Cc: bioconductor at stat.math.ethz.ch > Sent: Wednesday, March 23, 2011 6:14:50 PM > Subject: Re: [BioC] limma neqc > > Hi Ina: > > You will have to provide the "status" parameter to neqc then. This is a character vector which tells neqc which probes are regular one, which are negative control and which are positive controls. Regular probes should have a value of "regular" (same with value of the "regular" parameter) and negative controls should have a value of "negative". For positive controls, you can use their original names like "Housekeeping", "Biotin" etc, or you can just use a value of "positive" for all of them. The output of neqc will be a matrix which only includes regular probes. > > Please let me know if you have any further questions. > > Cheers, > Wei > > > On Mar 24, 2011, at 1:15 AM, Ina Hoeschele wrote: > >> Hi Wei, >> thanks much for your quick response. However, I should have mentioned that I am not reading the data using read.ilmn, I am reading and pre-processing the bead-level data via beadarray (doing QC etc.), then summarizing the data and then calling neqc. >> Ina >> >> ----- Original Message ----- >> From: "Wei Shi" <shi at="" wehi.edu.au=""> >> To: "Ina Hoeschele" <inah at="" vbi.vt.edu=""> >> Cc: bioconductor at stat.math.ethz.ch >> Sent: Tuesday, March 22, 2011 6:08:34 PM >> Subject: Re: [BioC] limma neqc >> >> Dear Ina: >> >> You do not need to specify positive control probes. They were automatically read in when you use read.ilmn function to read in both regular probe data and control probe data. You can use the code below to see the numbers of regular probes, negative control probes and positive control probes: >> >> x <- read.ilmn(files="probe profile.txt",ctrlfiles="control probe profile.txt", other.columns="Detection") # you will need to changes file names to the names of your files >> table(x$genes$Status) >> >> The neqc function will then perform a quantile normalization for data included in x using all the probes (including regular probes, negative controls and positive controls). The two parameters negctrl and regular tell neqc what the identifiers of regular probes and negative control probes are. The negative control probes play a particular role in the background correction (normexp background correction aided by negative controls). By default, regular probes have an identifier of "regular". Negative control probes by default have an identifier of "negative" which was the identifier used by almost all the BeadChip data we have ever seen. The neqc function should work in most cases with its default setting. >> >> For more details, please refer to the limma user guide for analyzing Illumina BeadChip data (section 11.7). >> >> Hope this helps. >> >> Cheers, >> Wei >> >> On Mar 23, 2011, at 6:40 AM, Ina Hoeschele wrote: >> >>> Hi Gordon et al., >>> I am confused about something - I am using limma neqc and according to a recent paper (Shi, Oshlack, Smyth 2010) the quantile normalization should be performed using both positive and negative control probes. But in the arguments for neqc one can specify the identifier for negative control probes and "regular" probes. So I set negctrl="negative" and regular="Gene". But this does not allow me to specify positive control probes? >>> Many thanks, Ina >>> >>> _______________________________________________ >>> Bioconductor mailing list >>> Bioconductor at r-project.org >>> https://stat.ethz.ch/mailman/listinfo/bioconductor >>> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor >> >> >> ______________________________________________________________________ >> The information in this email is confidential and intended solely for the addressee. >> You must not disclose, forward, print or use it without the permission of the sender. >> ______________________________________________________________________ > > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:14}}
ADD REPLYlink written 8.4 years ago by Wei Shi3.1k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 334 users visited in the last hour