IBG Wiki

User Tools

Site Tools

Trace:
keller_and_evans_lab:gscan_db_ga_p

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
gscan_db_ga_p [2017/02/14 04:35] – /* Phenotypes */ scottkeller_and_evans_lab:gscan_db_ga_p [2019/10/31 18:28] (current) – ↷ Links adapted because of a move operation 66.249.87.23
Line 576: Line 576:
  
 ====== ARIC ====== ====== ARIC ======
- 
-(Hannah/Joyce to update this section following Framingham as a guide) 
- 
-===== ID Mapping ===== 
  
  
Line 864: Line 860:
 ### sort ARIC.AFR.covariates.ped | uniq > ARIC.AFR.covariates.ped ### sort ARIC.AFR.covariates.ped | uniq > ARIC.AFR.covariates.ped
  
- 
- 
- 
-===== Genotypes ===== 
  
  
Line 878: Line 870:
 ===== Phenotypes ===== ===== Phenotypes =====
  
-Description of phenotypes can be found here: {{file_mesa_phenotypes_-_final.pdf}}+Description of phenotypes can be found here: {{:file_mesa_phenotypes_-_final.pdf}}
  
  
 ====== eMERGE ====== ====== eMERGE ======
  
-(Hannah/Joyce to update following Framingham as a guide) 
  
  
 ===== Phenotypes ===== ===== Phenotypes =====
  
-Description of phenotypes can be found here: {{file_emerge.pdf}}+ 
 + 
 +options(stringsAsFactors=F) 
 + 
 + 
 +### eMERGE is broken into different consent classes. We can conduct analyses on hmb, hmb-gso-nic, and  
 + 
 + 
 +emerge.hmb <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c1.HMB/PhenotypeFiles/phs000360.v3.pht003255.v2.p1.c1.MergedSet_Subject_Phenotypes.HMB.txt.gz", header=TRUE, sep="\t", stringsAsFactors=F) 
 +emerge.hmb.genos <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c1.HMB/GenotypeFiles/matrix/c1.HMB/eMerge_660_11212012_c1.fam", header=FALSE, sep="\t", stringsAsFactors=F) 
 + 
 + 
 +emerge.hmb.gso.nic <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c3.HM-B-GSO-NIC/PhenotypeFiles/phs000360.v3.pht003255.v2.p1.c3.MergedSet_Subject_Phenotypes.HM-B-GSO-NIC.txt.gz", header=TRUE, sep="\t", stringsAsFactors=F) 
 +emerge.hmb.gso.nic.genos <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c3.HM-B-GSO-NIC/GenotypeFiles/matrix/c3.HM-B-GSO-NIC/eMerge_660_11212012_c3.fam", header=FALSE, sep="\t", stringsAsFactors=F) 
 + 
 + 
 +emerge.hmb.gso <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c4.HMB-GSO/PhenotypeFiles/phs000360.v3.pht003255.v2.p1.c4.MergedSet_Subject_Phenotypes.HMB-GSO.txt.gz", header=TRUE, sep="\t", stringsAsFactors=F) 
 +emerge.hmb.gso.genos <- read.table("/work/KellerLab/dbGaP/eMERGE-MergedSet/PhenoGenotypeFiles/RootStudyConsentSet_phs000360.eMERGE_MergedSet.v3.p1.c4.HMB-GSO/GenotypeFiles/matrix/c4.HMB-GSO/eMerge_660_11212012_c4.fam", header=FALSE, sep="\t", stringsAsFactors=F) 
 + 
 + 
 +### Merge all files above according to SUBJID, which is used in the 
 +### genotype files. 
 + 
 +emerge <- merge(emerge.hmb, emerge.hmb.gso, all=T) 
 +emerge <- merge(d, emerge.hmb.gso.nic, all=T) 
 + 
 +### SMOKING INITIATION 
 +### 
 +### The eMERGE variable name is SMOKING_STATUS 
 +###      C65108 = never smoker 
 +###      C67147 = current smoker 
 +###      C67148 = past smoker 
 +###      C67151 = Unknown if ever smoked 
 +### 
 +### Descriptives: 
 +### 
 +### table(emerge$SMOKING_STATUS) 
 +### 
 +### C65108 C67147 C67148 C67151  
 +###   2217   1736   3457   9635  
 + 
 +si <- emerge$SMOKING_STATUS 
 +si[si == "C67147" | si == "C67148"] <- 2 
 +si[si == "C65108"] <- 1 
 +si[si != 1 & si != 2] <- NA 
 + 
 +### SMOKING Cessation 
 +### 
 +### Current == 2 & Former == 1 in GSCAN. This is already the case for these data. 
 + 
 +sc <- emerge$SMOKING_STATUS 
 +sc[sc == "C67147"] <- 2 
 +sc[sc == "C67148"] <- 1 
 +sc[sc != 1 & sc != 2] <- NA 
 + 
 + 
 +### eMERGE age variable is tricky because there is no obvious age at 
 +### assessment. We will use their "DECADE_BIRTH" as a terrible 
 +### approximation. 
 +### 1=1900-1919; 2=1920-1929, 3=1930-1939; 4=1940-1949; 5=1950-1959; 6=Unknown 
 +### 
 +### Descriptives: 
 +### 
 +### table(emerge$DECADE_BIRTH) 
 +### 
 +###   .    1    2    3    4    5    6    7    8    9   99  
 +###   6  612 2667 3533 4439 3127 1291  761  490   10  109  
 +birthyear <- emerge$DECADE_BIRTH 
 +birthyear[birthyear == "99"] <- NA 
 +birthyear[birthyear == "."] <- NA 
 + 
 + 
 +### SEX 
 +sex <- emerge$SEX 
 +sex[sex == "C46109"] <- 1 
 +sex[sex == "C46110"] <- 2 
 + 
 + 
 +### Scott decided not to correct for additional case-control variables 
 +### given what appears to be a highly complex sample and uncertainty 
 +### about the best course of action to account for disease status in 
 +### conducting smoking analyses. 
 + 
 +phenotypes <- data.frame(fid = emerge$SUBJID, 
 +                         iid = emerge$SUBJID, 
 +                         patid = "x", 
 +                         matid = "x", 
 +                         sex = sex, 
 +                         si = si, 
 +                         sc = sc) 
 + 
 +phenotypes[is.na(phenotypes)] <- "x" 
 + 
 +write.table(phenotypes, 
 +            "/work/KellerLab/vrieze/GSCAN/GWAS/summary_stats_generated_internally/eMERGE/GSCAN_eMERGE_phenotypes.ped", 
 +            row.names=F, 
 +            quote = F, 
 +            sep="\t") 
 + 
 + 
 +covariates  <- data.frame(fid = emerge$SUBJID, 
 +                          iid = emerge$SUBJID, 
 +                          patid = "x", 
 +                          matid = "x", 
 +                          sex = sex, 
 +                          birthyear = birthyear)                          
 + 
 +covariates[is.na(covariates)] <- "x" 
 + 
 +write.table(covariates, 
 +            "/work/KellerLab/vrieze/GSCAN/GWAS/summary_stats_generated_internally/eMERGE/GSCAN_eMERGE_covariates.ped", 
 +            row.names=F, 
 +            quote = F, 
 +            sep="\t") 
  
  
Line 1040: Line 1145:
  
  
-beagess_data <- read.table("/work/KellerLab/GSCAN/dbGaP/Barretts/PhenoGenotypeFiles/RootStudyConsentSet_phs000869.BEAGESS.v1.p1.c1.GRU-MDS/PhenotypeFiles/phs000869.v1.pht004610.v1.p1.c1.BEAGESS_Subject_Phenotypes.GRU-MDS.txt.gz", header=TRUE, sep="\t", stringsAsFactors=F) +options(stringsAsFactors=F)
-geno_data <- read.table("/work/KellerLab/GSCAN/dbGaP/Barretts/PhenoGenotypeFiles/RootStudyConsentSet_phs000869.BEAGESS.v1.p1.c1.GRU-MDS/GenotypeFiles/phg000580.v1.NCI_BEAGESS.genotype-calls-matrixfmt.c1.GRU-MDS.update/BEAGESS_dbGaP_29Jan2015.fam", header=FALSE, stringsAsFactors=F) +
-geno_data <- geno_data[-c(6)] +
-geno_data$V3[geno_data$V3 =="0"] = "x" +
-geno_data$V4[geno_data$V4 =="0"] = "x"+
  
-### IMPORTANT: All columns in original phenotype file are shifted one column to the the right, so the label of the column does not match the data in that column!!! +beagess <- read.table("/work/KellerLab/GSCAN/dbGaP/Barretts/PhenoGenotypeFiles/RootStudyConsentSet_phs000869.BEAGESS.v1.p1.c1.GRU-MDS/PhenotypeFiles/phs000869.v1.pht004610.v1.p1.c1.BEAGESS_Subject_Phenotypes.GRU-MDS.txt.gz", header=TRUE, sep="\t", stringsAsFactors=F)
-si <- beagess_data$BMI +
-sc <- beagess_data$BMI +
-age <- beagess_data$sex+
  
-beagess_data <- cbind(beagess_data, beagess_data[,6]) +genos <- read.table("/work/KellerLab/GSCAN/dbGaP/Barretts/PhenoGenotypeFiles/RootStudyConsentSet_phs000869.BEAGESS.v1.p1.c1.GRU-MDS/GenotypeFiles/phg000580.v1.NCI_BEAGESS.genotype-calls-matrixfmt.c1.GRU-MDS.update/BEAGESS_dbGaP_29Jan2015.fam", header=FALSE, stringsAsFactors=F)
-colnames(beagess_data)[13] <- "sc"+
  
  
-### SMOKING INITIATION +### The file reads into R incorrectly because of a weird trailing tab 
 +### character in the data file, so use the below code to shift column 
 +### names to the correct column. 
 +names(beagess)[1] <- "XXX" 
 +beagess$dbGaP_Subject_ID <- row.names(beagess) 
 +beagess$assocEABEvsCO <- NULL 
 +names(beagess) <- c(names(beagess)[2:11], "assocEABEvsCO", "dbGaP_Subject_ID") 
 + 
 + 
 +### SMOKING INITIATION
 ### ###
 ### BEAGESS variable name is "cig_smk_status" and is listed under the column "BMI". ### BEAGESS variable name is "cig_smk_status" and is listed under the column "BMI".
Line 1062: Line 1168:
 ###       -99 = not consented ###       -99 = not consented
 ###       -9 = Missing ###       -9 = Missing
-###        0 = Never +###        0 = Never
 ###        1 = Former ###        1 = Former
 ###        2 = Current ###        2 = Current
-###     +### 
-### Descriptives: +### Descriptives:
 ### ###
 ### > table(beagess_data$BMI) ### > table(beagess_data$BMI)
-###  +### 
-### -99   -9    0    1    2  +### -99   -9    0    1    2 
-### 494 2051 1515 2288  575 +### 494 2051 1515 2288  575
 ### ###
 ### > summary(beagess_data$BMI) ### > summary(beagess_data$BMI)
-### Min. 1st Qu.  Median    Mean 3rd Qu.    Max.  +### Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-### -99.000  -9.000   0.000  -9.234   1.000   2.000 +### -99.000  -9.000   0.000  -9.234   1.000   2.000
  
-beagess_data$BMI[beagess_data$BMI == -99] = "x" +si <- beagess$cig_smk_status 
-beagess_data$BMI[beagess_data$BMI == -9] = "x" +si[si == 1 | si == 2] <- 2 
-beagess_data$BMI[beagess_data$BMI == 1] = 2 +si[si == 0] <- 1 
-beagess_data$BMI[beagess_data$BMI == 0] = 1+si[si != 1 & si != 2] <- NA
  
-### SMOKING Cessation +### SMOKING Cessation
 ### ###
 ### BEAGESS variable name is "cig_smk_status" and is listed under the column "BMI". ### BEAGESS variable name is "cig_smk_status" and is listed under the column "BMI".
Line 1089: Line 1195:
 ###       -99 = not consented ###       -99 = not consented
 ###       -9 = Missing ###       -9 = Missing
-###        0 = Never +###        0 = Never
 ###        1 = Former ###        1 = Former
 ###        2 = Current ###        2 = Current
-###  
-### Descriptives:  
 ### ###
-### > table(beagess_data$BMI) +### Descriptives: 
-###  +### table(beagess$cig_smk_status) 
-### -99   -9    0    1    2  +###  -99   -9    0    1    2  
-### 494 2051 1515 2288  575 +###  494 2051 1515 2288  575 
 ### ###
-### > summary(beagess_data$BMI) +### Current == 2 & Former == 1 in GSCAN. This is already the case for these data.
-### Min. 1st Qu.  Median    Mean 3rd Qu.    Max.  +
-### -99.000  -9.000   0.000  -9.234   1.000   2.000 +
  
-beagess_data$sc[beagess_data$sc == -99] = "x" +sc <- beagess$cig_smk_status 
-beagess_data$sc[beagess_data$sc == -9] = "x" +sc[sc == 0 | sc == -9 | sc == -99] <- NA
-beagess_data$sc[beagess_data$sc == 0] = "x"+
  
-### AGE  +### BEAGESS variable name is "agegpcat"
-### +
-### BEAGESS variable name is "agegpcat" and is listed under the column "sex".+
 ###   "Age in 5 year categories" ###   "Age in 5 year categories"
-###   Response options are integers 1 through 14 or -9. +###   Response options are integers 1 through 14 or -9.
 ###       -9 = Missing ###       -9 = Missing
 ###       1 = 15-29 years of age ###       1 = 15-29 years of age
Line 1129: Line 1228:
 ###       14 = 90-100 years of age ###       14 = 90-100 years of age
 ### ###
-### Descriptives: +### Descriptives:
 ### ###
-### > table(beagess_data$sex) +### > table(beagess_data$agegpcat)
-###  +
-### -9    1    2    3    4    5    6    7    8    9   10   11   12   13   14  +
-### 27   19   48  112  190  378  657  958 1134 1238 1018  794  246   87   17 +
 ### ###
-### > summary(beagess_data$BMI) +### -9    1    2    3    4    5    6    7    8    9   10   11   12   13   14 
-### Min. 1st Qu.  Median    Mean 3rd Qu.    Max.  +### 27   19   48  112  190  378  657  958 1134 1238 1018  794  246   87   17 
-### -9.000   7.000   8.000   8.227  10.000  14.000 +### 
 +### This is fine, but change the -9 to NA 
 +beagess$agegpcat[beagess$agegpcat == -9] <- NA
  
-beagess_data$sex[beagess_data$sex== 1] = 22 
-beagess_data$sex[beagess_data$sex== 2] = 32 
-beagess_data$sex[beagess_data$sex== 3] = 37 
-beagess_data$sex[beagess_data$sex== 4] = 42 
-beagess_data$sex[beagess_data$sex== 5] = 47 
-beagess_data$sex[beagess_data$sex== 6] = 52 
-beagess_data$sex[beagess_data$sex== 7] = 57 
-beagess_data$sex[beagess_data$sex== 8] = 62 
-beagess_data$sex[beagess_data$sex== 9] = 67 
-beagess_data$sex[beagess_data$sex== 10] = 72 
-beagess_data$sex[beagess_data$sex== 11] = 77 
-beagess_data$sex[beagess_data$sex== 12] = 82 
-beagess_data$sex[beagess_data$sex== 13] = 87 
-beagess_data$sex[beagess_data$sex== 14] = 95 
-beagess_data$sex[beagess_data$sex== -9] = "x" 
-beagess_data$sex <- as.numeric(beagess_data$sex) 
  
-### rename SUBJID column to dbGaP_Subject_ID in genotype file so it matches phenotype file where SUBJID is misslabelled as dbGaP_Subject_ID 
-colnames(geno_data) <- c("famid", "dbGaP_Subject_ID", "patid", "matid", "sex") 
  
-### merge geno and pheno files +### looks like BEAGESS uses their own internal "SUBJECT_ID" in the 
-phen <- merge(geno_data,beagess_data, by="dbGaP_Subject_ID", all =TRUE) +### genotype files, so we'll use that in our pedigree files 
-phen <- phen[,c(1:5,7,10,17)]+phenotypes <- data.frame(fid = beagess$SUBJECT_ID, 
 +                         iid = beagess$SUBJECT_ID, 
 +                         patid = "x", 
 +                         matid = "x", 
 +                         sex = beagess$sex, 
 +                         si = si, 
 +                         sc = sc)
  
-phen <- phen[,c(2,1,3,4,5,6,7,8)] +phenotypes[is.na(phenotypes)] <- "x"
-colnames(phen) <- c("famid", "dbGaP_Subject_ID", "patid", "matid", "sex", "age","si","sc")+
  
 +write.table(phenotypes,
 +            "/work/KellerLab/vrieze/GSCAN/GWAS/summary_stats_generated_internally/BEAGESS/GSCAN_BEAGESS_phenotypes.ped",
 +            row.names=F,
 +            quote = F,
 +            sep="\t")
  
-phenotypes <- data.frame(famid=phen$famid, dbGaP_Subject_ID=phen$dbGaP_Subject_ID, patid=phen$patid, matid=phen$matid, sex=phen$sex, si=phen$si, sc=phen$sc, currentformersmoker=phen$sc, age=phen$age, age2=phen$age^2) 
-colnames(phenotypes) <- c("famid", "dbGaP_Subject_ID", "patid", "matid", "sex", "si", "sc", "currentformersmoker", "age", "age2") 
  
-pcs <- read.table("/rc_scratch/hayo0753/BEAGESS/BEAGESS_pcs_and_ancestrys.txt",head=TRUE, stringsAsFactors=F, sep=" ") +covariates  <- data.frame(fid = beagess$SUBJECT_ID, 
-colnames(pcs) [1] <- "dbGaP_Subject_ID"+                          iid = beagess$SUBJECT_ID, 
 +                          patid = "x", 
 +                          matid = "x", 
 +                          sex = beagess$sex, 
 +                          age = beagess$age, 
 +                          age2 = beagess$age^2, 
 +                          Barretts.esophagus.case_v_control = beagess$assocBEvsCO, 
 +                          Esophageal.adenocarcinoma.case_v_control = beagess$assocEAvsCO)
  
-gscan.phenotypes <- merge(phenotypes, pcs, by="dbGaP_Subject_ID", all=TRUE) +covariates[is.na(covariates)] <- "x"
-gscan.phenotypes[is.na(gscan.phenotypes)] <- "x"+
  
-### EUROPEANS - entire sample is EUR +write.table(covariates, 
-phenotypes.EUR.ped <- subset(gscan.phenotypes, ancestry == "EUR", select=c("famid","dbGaP_Subject_ID","patid","matid", "sex", "si", "sc")) +                        "/work/KellerLab/vrieze/GSCAN/GWAS/summary_stats_generated_internally/BEAGESS/GSCAN_BEAGESS_covariates.ped", 
-write.table(phenotypes.EUR.ped, file="BEAGESS.EUR.phenotypes.ped", quote=F, col.names=T, row.names=F, sep="\t") +                        row.names=F, 
- +                        quote = F, 
-covariates.EUR.ped <- subset(gscan.phenotypes, ancestry == "EUR", select=c("famid","dbGaP_Subject_ID","patid","matid", "sex", "age", "age2",  "currentformersmoker", "PC1", "PC2", "PC3", "PC4", "PC5", "PC6", "PC7", "PC8", "PC9", "PC10")) +                        sep="\t") ### save table                 
-write.table(covariates.EUR.ped, file="BEAGESS.EUR.covariates.ped", quote=F, col.names=T, row.names=F, sep="\t")+
  
  
keller_and_evans_lab/gscan_db_ga_p.1487046907.txt.gz · Last modified: by scott