library(RODBC)
library(data.table)
options("width"=400)
dbCohort="VA_MVP011_CDW.[KDI\\mcmahonb.mvp011w].Cc2022"  #change also connection string ch1
dbProj="VA_MVP011_CDW"
ndbCohort="MVP_CDWWork.BenMcMahonTemp.[KDI\\mcmahonb.mvp011w].Cc2022"  #uses temporaty table instead of linking across servers.

ch2 <- odbcDriverConnect('Driver={ODBC Driver 13 for SQL Server};Server=vawnornl.va.ornlkdi.org;Trusted_Connection=yes;ServerSPN=MSSQLSvc/vawnornl.va.ornlkdi.org') 
ch <- odbcDriverConnect('Driver={ODBC Driver 13 for SQL Server};Server=vawnornl.va.ornlkdi.org;Database=CDWWork;Trusted_Connection=yes;ServerSPN=MSSQLSvc/vawnornl.va.ornlkdi.org') 
ch1=odbcDriverConnect(paste(c('Driver={ODBC Driver 13 for SQL Server};Server=vawnornl.va.ornlkdi.org;Database=',dbProj,';Trusted_Connection=yes;ServerSPN=MSSQLSvc/vawnornl.va.ornlkdi.org'),collapse=''))

sqlQuery(ch2,query="select top 10 * from MVP_CDWWork.CDWWork.Outpat.VProcedure")
sqlQuery(ch2,query=paste(c("select top 10 * from ",dbProj,".SPatient.SPatient"),collapse=""))

sqlQuery(ch2,query=paste(c("select top 10 * from A.PatientICN, A.BirthDateTime, B.Sta3n from ",dbProj,".SPatient.SPatient A INNER JOIN MVP_CDWWork.CDWWork.Outpat.VProcedure B on A.PatientSID=B.PatientSID"),collapse=""))

tbls=sqlTables(ch)
table(tbls$TABLE_SCHEM)
tbls[tbls$TABLE_SCHEM == "Dim",]$TABLE_NAME

Mode=function(x){ux=unique(x);ux[which.max(tabulate(match(x,ux)))]}

dob=sqlQuery(ch2,query="select PatientICN, PatientSID,MOB,YOB from VA_MVP011_CDW.SPatient.SPatient",stringsAsFactors=F) 
dob=data.table(dob)               # provides 24,860,108 unique patients, 18 January, 2022
dob=dob[is.na(YOB)==F,]           # provides 24,834,168 unique patients
dob=dob[YOB>=1901 & YOB <=2003,]  #leaves 24,558,159 unique patients
pts=unique(dob,by="PatientICN")
pts$day=runif(dim(pts)[1],0,28)
pts$YoB=pts$YOB-2000+(pts$MOB+pts$day/30.4375)/12

setkey(pts,YoB)
pts$RN=seq.int(nrow(pts))
pts$C=ceiling(pts$RN/245600)-1

birth=data.table(pts$PatientICN,pts$YoB,pts$C)
colnames(birth)=c("PatientICN","YoB","C")
setkey(dob,PatientICN)
setkey(birth,PatientICN)
cohorts=merge(dob,birth,all.x=T,by="PatientICN")
cohorts=cohorts[is.na(PatientICN)==F,]
cohorts$MOB=NULL;cohorts$YOB=NULL

saveRDS(cohorts,file="Cohort.2020.RDS")
sqlSave(ch1,C,tablename="Ci2020",rownames=F)  # took 5hrs 17m

#minimal demographic info n (written to ../tmp, to be modified twice below, finally called DemoU.xx.RDS)
for (cohort in seq(0,100)){
q=paste(c("select A.PatientICN, A.BirthDateTime, A.DeceasedFlag, A.DeathDateTime, A.Gender, A.MaritalStatus from ",dbProj,".SPatient.SPatient A INNER JOIN ",ndbCohort," B ON A.PatientICN = B.PatientICN WHERE B.C =",cohort),collapse="")
C=sqlQuery(ch2,query=q,stringsAsFactors=F)
C=data.table(C)
d=unique(C, by="PatientICN")
C=data.table(d)
C[,DeathDateTime:= as.numeric(julian(DeathDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(C,PatientICN)
print(c(dim(C),cohort))
saveRDS(C,file=paste(c("Demo.",cohort,".RDS"),collapse=""))
}

# race from PatSub.PatientRace table
for (cohort in seq(0,100)){
q=paste(c("select A.PatientRaceSID, C.PatientICN, A.PatientSID, A.RaceSID, A.Race, A.CollectionMethod from MVP_CDWWork.CDWWork.PatSub.PatientRace A INNER JOIN ",ndbCohort," C ON A.PatientSID = C.PatientSID WHERE C.C =",cohort),collapse="")
R=sqlQuery(ch2,query=q,stringsAsFactors=F)
R=data.table(R)
Rs=R[CollectionMethod=="SELF IDENTIFICATION",Mode(Race),by="PatientICN"]
Ro=R[is.na(match(PatientICN,Rs$PatientICN)),Mode(Race),by="PatientICN"]
Rk=rbind(Rs,Ro)
colnames(Rk)[2]="Race"
setkey(Rk,PatientICN)
# ethnicity from PatSub.Ethnicity table
q=paste(c("select A.PatientEthnicitySID, C.PatientICN, A.PatientSID, A.Ethnicity, A.CollectionMethod from MVP_CDWWork.CDWWork.PatSub.PatientEthnicity A INNER JOIN ",ndbCohort," C ON A.PatientSID = C.PatientSID WHERE C.C =",cohort),collapse="")
E=sqlQuery(ch2,query=q,stringsAsFactors=F)
E=data.table(E)
Es=E[CollectionMethod=="SELF IDENTIFICATION",Mode(Ethnicity),by="PatientICN"]
Eo=E[is.na(match(PatientICN,Es$PatientICN)),Mode(Ethnicity),by="PatientICN"]
Ek=rbind(Es,Eo)
colnames(Ek)[2]="Ethnicity"
setkey(Ek,PatientICN)
demo=readRDS(,file=paste(c("Demo.",cohort,".RDS"),collapse=""))
demo2=Rk[demo]
demo3=Ek[demo2]
demo4=demo3[,c(1,4,5,6,7,8,3,2)]
print(c(dim(demo4),cohort))
saveRDS(demo4,file=paste(c("DemoRE.",cohort,".RDS"),collapse=""))    #(written) to ../tmp/, be modified below with healthcare utilization data, saved to DemoU.xx.RDS
}

# ontologies and labels
dimP9=sqlQuery(ch,query="select A.ICD9ProcedureCode, B.ICD9ProcedureDescription from Dim.ICD9Procedure A INNER JOIN Dim.ICD9ProcedureDescriptionVersion B ON A.ICD9ProcedureSID=B.ICD9ProcedureSID",stringsAsFactors=F)
dimP9=data.table(dimP9)
P9=unique(dimP9,by="ICD9ProcedureCode")
saveRDS(P9,file="Dim.ICD9Proc.RDS")

P10=sqlQuery(ch,query="select  A.ICD10ProcedureCode, B.ICD10ProcedureDescription from Dim.ICD10Procedure A INNER JOIN Dim.ICD10ProcedureDescriptionVersion B on A.ICD10ProcedureSID = B.ICD10ProcedureSID",stringsAsFactors=F)
P10 = data.table(P10)
Pr10=unique(P10, by="ICD10ProcedureCode")
P10 = data.table(P10)
setkey(Pr10,ICD10ProcedureCode)
saveRDS(P10,file="Dim.ICD10Proc.RDS")

CPT=sqlQuery(ch,query="select A.CPTCode, A.CPTName, A.CPTDescription, B.CPTCategory, B.MajorCPTCategory from Dim.CPT A LEFT JOIN Dim.CPTCategory B on A.CPTCategorySID=B.CPTCategorySID",stringsAsFactors=F)
CPT=data.table(CPT)
CPT=unique(CPT,by="CPTCode")
saveRDS(CPT,file="Dim.CPT.RDS")

Drug=sqlQuery(ch,query="select DrugClassCode, DrugClassification from Dim.DrugClass",stringsAsFactors=F)
Drug=data.table(Drug)
Dr=unique(Drug,by="DrugClassCode")
setkey(Dr,DrugClassCode)
saveRDS(Dr,file="Dim.Drugs.RDS")

DimMed18=sqlQuery(ch1,query="select * from Dim.DimMed_v18_2",stringsAsFactors=F) 
DimMed18=data.table(DimMed18)
setkey(DimMed18,LocalDrugSID)
saveRDS(DimMed18,file="DimMed18.RDS")

DimLab18=sqlQuery(ch1,query="select * from Dim.DimLab_v18_2",stringsAsFactors=F) 
DimLab18=data.table(DimLab18)
setkey(DimLab18,DimLabSID)
saveRDS(DimLab18,file="DimLab18.RDS")

Di9=sqlQuery(ch,query="select  B.ICDCode, A.ICDDescription from Dim.ICD A INNER JOIN Dim.ICD B on A.ICDSID = B.ICDSID",stringsAsFactors=F)
Di9=data.table(Di9)
D9=unique(Di9,by="ICDCode")
setkey(D9,ICDCode)
saveRDS(D9,file="Dim.ICD9.RDS")

Di10=sqlQuery(ch,query="select  B.ICD10Code, A.ICD10Description from Dim.ICD10DescriptionVersion A INNER JOIN Dim.ICD10 B on A.ICD10SID = B.ICD10SID",stringsAsFactors=F)
Di10=data.table(Di10)
D10=unique(Di10,by="ICD10Code")
setkey(D10,ICD10Code)
saveRDS(D10,file="Dim.ICD10.RDS")

V=sqlQuery(ch,query="select VitalType, VitalTypeAbbreviation from Dim.VitalType ",stringsAsFactors=F)
V=data.table(V)
V=unique(V,by="VitalType")
setkey(V,VitalType)
saveRDS(V,file="Dim.Vital.RDS")

LOINC=sqlQuery(ch,query="SELECT LOINC,Component,Property,TimeAspect,ScaleType,MethodType FROM Dim.LOINC",stringsAsFactors=F)
LOINC=data.table(LOINC)
LOINC=unique(LOINC,by="LOINC")
setkey(LOINC,LOINC)
saveRDS(LOINC,file="Dim.LOINC.RDS")

#  Vital signs
for (cohort in seq(0,99)){
q=paste(c("select C.PatientICN, A.VitalSignTakenDateTime, A.Systolic, A.Diastolic from MVP_CDWWork.CDWWork.Vital.VitalSign A INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.Systolic > 10 AND C.C = ",cohort),collapse="")
V=sqlQuery(ch2,query=q,stringsAsFactors=F)
V=data.table(V)
V[,VitalSignTakenDateTime:= as.numeric(julian(VitalSignTakenDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(V,PatientICN)
saveRDS(V,file=paste(c("Vital.bp.",cohort,".RDS"),collapse=""))
}

for (cohort in seq(0,99)){
q=paste(c("select C.PatientICN, A.VitalSignTakenDateTime, A.VitalResult, A.VitalResultNumeric, B.VitalTypeAbbreviation from MVP_CDWWork.CDWWork.Vital.VitalSign A INNER JOIN MVP_CDWWork.CDWWork.Dim.VitalType B on A.VitalTypeSID=B.VitalTypeSID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.VitalResultNumeric > 0.1 AND C.C = ",cohort),collapse="")
V=sqlQuery(ch2,query=q,stringsAsFactors=F)
V=data.table(V)
V[,VitalSignTakenDateTime:= as.numeric(julian(VitalSignTakenDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(V,PatientICN)
saveRDS(V,file=paste(c("Vital.other.",cohort,".RDS"),collapse=""))
}

#inpatient dx
for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.InpatientSID, A.Sta3n, A.DischargeDateTime, B.ICDCode, B.ICDDescription from MVP_CDWWork.CDWWork.Inpat.InpatientDiagnosis A INNER JOIN MVP_CDWWork.CDWWork.Dim.ICD B ON A.ICD9SID=B.ICDSID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.ICD9SID != 0 AND C.C =",cohort),collapse="")
icd9=sqlQuery(ch2,query=q,stringsAsFactors=F)
icd9=data.table(icd9)
icd9[,DischargeDateTime:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(icd9,PatientICN)
saveRDS(icd9,file=paste(c("iDx9.",cohort,".RDS"),collapse=""))
print(c(dim(icd9),cohort))
}

for (cohort in seq(1,100)){   #cohort 0 has no icd10 inpatient Dx codes
q=paste(c("select C.PatientICN, A.InpatientSID, A.Sta3n, A.DischargeDateTime, B.ICD10Code from MVP_CDWWork.CDWWork.Inpat.InpatientDiagnosis A INNER JOIN MVP_CDWWork.CDWWork.Dim.ICD10 B ON A.ICD10SID=B.ICD10SID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.ICD10SID != 0 AND C.C =",cohort),collapse="")
#q=paste(c("select C.PatientICN, A.InpatientSID, A.Sta3n, A.DischargeDateTime, B.ICD10Code from MVP_CDWWork.CDWWork.Inpat.InpatientDiagnosis A INNER JOIN MVP_CDWWork.CDWWork.Dim.ICD10 B ON A.ICD10SID=B.ICD10SID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.ICD10SID != 0 AND C.C =",cohort),collapse="")
icd10=sqlQuery(ch2,query=q,stringsAsFactors=F)
icd10=data.table(icd10)
icd10[,DischargeDateTime:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(icd10,PatientICN)
saveRDS(icd10,file=paste(c("iDx10.",cohort,".RDS"),collapse=""))
print(c(dim(icd10),cohort))
}
cohort=40
for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.InpatientSID, A.Sta3n, A.AdmitDateTime, A.DischargeDateTime, B.ICDCode, B.ICDDescription from MVP_CDWWork.CDWWork.Inpat.InpatientDischargeDiagnosis A INNER JOIN MVP_CDWWork.CDWWork.Dim.ICD B ON A.ICD9SID=B.ICDSID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.ICD9SID != 0 AND C.C =",cohort),collapse="")
icd9=sqlQuery(ch2,query=q,stringsAsFactors=F)
icd9=data.table(icd9)
icd9[,AdmitDateTime:=as.numeric(julian(AdmitDateTime,origin=as.Date("2000-01-01")))/365.25]
icd9[,duration:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))-AdmitDateTime*365.25]
icd9[,DischargeDateTime:=NULL]
setkey(icd9,PatientICN)
saveRDS(icd9,file=paste(c("idDx9.",cohort,".RDS"),collapse=""))
print(c(dim(icd9),cohort))
}

for (cohort in seq(1,100)){   #cohort 0 has no icd10 inpatient Dx codes
q=paste(c("select C.PatientICN, A.InpatientSID, A.Sta3n, A.AdmitDateTime, A.DischargeDateTime, B.ICD10Code from MVP_CDWWork.CDWWork.Inpat.InpatientDischargeDiagnosis A INNER JOIN MVP_CDWWork.CDWWork.Dim.ICD10 B ON A.ICD10SID=B.ICD10SID INNER JOIN ",ndbCohort," C on A.PatientSID=C.PatientSID WHERE A.ICD10SID != 0 AND C.C =",cohort),collapse="")
icd10=sqlQuery(ch2,query=q,stringsAsFactors=F)
icd10=data.table(icd10)
icd10[,AdmitDateTime:=as.numeric(julian(AdmitDateTime,origin=as.Date("2000-01-01")))/365.25]
icd10[,duration:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))-AdmitDateTime*365.25]
icd10[,DischargeDateTime:=NULL]
setkey(icd10,PatientICN)
saveRDS(icd10,file=paste(c("idDx10.",cohort,".RDS"),collapse=""))
print(c(dim(icd10),cohort))
}
#inpatient discharge types
for (cohort in seq(0,100)){
    q=paste(c("select C.PatientICN, A.InpatientSID, A.DischargeDateTime, A.DispositionType, B.Specialty, B.PTFCode, D.PlaceOfDisposition from MVP_CDWWork.CDWWork.Inpat.Inpatient A LEFT JOIN MVP_CDWWork.CDWWork.Dim.Specialty B on A.DischargeFromSpecialtySID=B.SpecialtySID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.PlaceOfDisposition D on A.PlaceOfDispositionSID=D.PlaceOfDispositionSID WHERE C.C = ",cohort),collapse="")
    iDtype=sqlQuery(ch2,query=q,stringsAsFactors=F)
    iDtype=data.table(iDtype)
    iDtype[,DischargeDateTime:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))/365.25]
    setkey(iDtype,PatientICN)
    saveRDS(iDtype,file=paste(c("iDtype.",cohort,".RDS"),collapse=""))
    print(c(dim(iDtype),cohort))
}

#inpatient ICD 9/10 procedure  
for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.InpatientSID, A.sta3n, A.AdmitDateTime, A.DischargeDateTime, B.ICD9ProcedureCode from MVP_CDWWork.CDWWork.Inpat.InpatientICDProcedure A INNER Join MVP_CDWWork.CDWWork.Dim.ICD9Procedure B on A.ICD9ProcedureSID = B.ICD9ProcedureSID INNER Join ",ndbCohort," C on A.PatientSID = C.PatientSID  WHERE A.ICD9ProcedureSID != '0' AND C.C = ",cohort),collapse="")
Inpat.Proc9=sqlQuery(ch2,query=q,stringsAsFactors=F)
Inpat.Proc9=data.table(Inpat.Proc9)
Inpat.Proc9[,AdmitDateTime:=as.numeric(julian(AdmitDateTime,origin=as.Date("2000-01-01")))/365.25]
Inpat.Proc9[,duration:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))-AdmitDateTime*365.25]
Inpat.Proc9[,DischargeDateTime:=NULL]
setkey(Inpat.Proc9,PatientICN)
saveRDS(Inpat.Proc9,file=paste(c("iProc9.",cohort,".RDS"),collapse=""))
print(c(dim(Inpat.Proc9),cohort))
}

for (cohort in seq(0,99)){   #(no icd10 codes for cohort 0)
q=paste(c("select C.PatientICN, A.InpatientSID, A.sta3n, A.AdmitDateTime, A.DischargeDateTime, B.ICD10ProcedureCode from MVP_CDWWork.CDWWork.Inpat.InpatientICDProcedure A INNER Join MVP_CDWWork.CDWWork.Dim.ICD10Procedure B on A.ICD10ProcedureSID = B.ICD10ProcedureSID INNER Join ",ndbCohort," C on A.PatientSID = C.PatientSID  WHERE A.ICD10ProcedureSID != '0' AND C.C =",cohort),collapse="")
Inpat.Proc10=sqlQuery(ch2,query=q,stringsAsFactors=F)
Inpat.Proc10=data.table(Inpat.Proc10)
Inpat.Proc10[,AdmitDateTime:=as.numeric(julian(AdmitDateTime,origin=as.Date("2000-01-01")))/365.25]
Inpat.Proc10[,duration:=as.numeric(julian(DischargeDateTime,origin=as.Date("2000-01-01")))-AdmitDateTime*365.25]
Inpat.Proc10[,DischargeDateTime:=NULL]
setkey(Inpat.Proc10,PatientICN)
saveRDS(Inpat.Proc10,file=paste(c("iProc10.",cohort,".RDS"),collapse=""))
print(c(dim(Inpat.Proc10),cohort))
}

#inpatient CPT procedure
cohort=50
for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.InpatientSID, B.CPTCode, A.CPTProcedureDateTime from MVP_CDWWork.CDWWork.Inpat.InpatientCPTProcedure A INNER JOIN MVP_CDWWork.CDWWork.Dim.CPT B ON A.CPTSID = B.CPTSID INNER JOIN ",ndbCohort," C ON A.PatientSID=C.PatientSID WHERE C.C = ",cohort),collapse="")
iCPT=sqlQuery(ch2,query=q,stringsAsFactors=F)
iCPT=data.table(iCPT)
iCPT[,CPTProcedureDateTime:=as.numeric(julian(CPTProcedureDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(iCPT,PatientICN)
saveRDS(iCPT,file=paste(c("iCPT.",cohort,".RDS"),collapse=""))
print(c(dim(iCPT),cohort))
}

#pharmacy pickups
for (cohort in seq(1,99)){
for (type in c("A","CN","CV","HS")){
q=paste(c("select C.PatientICN, A.LocalDrugSID, A.LocalDrugNameWithDose, A.FillDateTime, A.ReleaseDateTime, A.ReturnedToStockDateTime, A.QtyNumeric, A.DaysSupply, B.VAClassification, B.NationalDrugNameWithDose, B.DrugNameWithoutDose from MVP_CDWWork.CDWWork.RxOut.RxOutpatFill A Inner JOIN MVP_CDWWork.CDWWork.Dim.LocalDrug B ON A.LocalDrugSID=B.LocalDrugSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE (B.VAClassification LIKE '",type,"%' AND C.C =",cohort,")"),collapse="")
Rx=sqlQuery(ch2,query=q,stringsAsFactors=F)
Rx=data.table(Rx)
Rx[,FillDateTime:=as.numeric(julian(FillDateTime,origin=as.Date("2000-01-01")))/365.25]
Rx[,ReleaseDateTime:=as.numeric(julian(ReleaseDateTime,origin=as.Date("2000-01-01")))/365.25]
Rx[,ReturnedToStockDateTime:=as.numeric(julian(ReturnedToStockDateTime,origin=as.Date("2000-01-01")))/365.25]
saveRDS(Rx,file=paste(c("Rx.",type,".",cohort,".RDS"),collapse=""))
print(c(dim(Rx),cohort,type))
}

q=paste(c("select C.PatientICN, A.LocalDrugSID, A.LocalDrugNameWithDose, A.FillDateTime, A.ReleaseDateTime, A.ReturnedToStockDateTime, A.QtyNumeric, A.DaysSupply, B.VAClassification, B.NationalDrugNameWithDose, B.DrugNameWithoutDose from MVP_CDWWork.CDWWork.RxOut.RxOutpatFill A Inner JOIN MVP_CDWWork.CDWWork.Dim.LocalDrug B ON A.LocalDrugSID=B.LocalDrugSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE (B.VAClassification NOT LIKE 'CN%' and B.VAClassification NOT LIKE 'CV%' and B.VAClassification NOT LIKE 'A%' and B.VAClassification NOT LIKE 'HS%' AND C.C =",cohort,")"),collapse="")
Rx=sqlQuery(ch2,query=q,stringsAsFactors=F)
Rx=data.table(Rx)
Rx[,FillDateTime:=as.numeric(julian(FillDateTime,origin=as.Date("2000-01-01")))/365.25]
Rx[,ReleaseDateTime:=as.numeric(julian(ReleaseDateTime,origin=as.Date("2000-01-01")))/365.25]
Rx[,ReturnedToStockDateTime:=as.numeric(julian(ReturnedToStockDateTime,origin=as.Date("2000-01-01")))/365.25]
saveRDS(Rx,file=paste(c("Rx.Oth.",cohort,".RDS"),collapse=""))
print(c(dim(Rx),cohort))
}

#BCMA meds 
for (cohort in seq(0,100)){
    q=paste(c("select C.PatientICN, A.ActionDateTime, B.VAClassification, D.DrugNameWithoutDose, B.LocalDrugNameWithDose, B.NationalDrugNameWithDose, A.DosesOrdered, A.DosesGiven, A.UnitOfAdministration from MVP_CDWWork.CDWWork.BCMA.BCMADispensedDrug A INNER JOIN MVP_CDWWork.CDWWork.Dim.LocalDrug B on A.LocalDrugSID=B.LocalDrugSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID left join MVP_CDWWork.CDWWork.Dim.DrugNameWithoutDose D on B.DrugNameWithoutDoseSID=D.DrugNameWithoutDoseSID WHERE C.C =",cohort),collapse="")
    iMeds=sqlQuery(ch2,query=q,stringsAsFactors=F)
    iMeds=data.table(iMeds)
    iMeds[,ActionDateTime:=as.numeric(julian(ActionDateTime,origin=as.Date("2000-01-01")))/365.25]
    saveRDS(iMeds,file=paste(c("BCMA.",cohort,".RDS"),collapse=""))
    print(c(dim(iMeds),cohort))
}

#FEE   
for (cohort in seq(0,100)){
    q=paste(c("SELECT C.PatientICN, A.Sta3n, B.InitialTreatmentDateTime, D.FeePurposeOfVisit, E.CPTCode FROM ",ndbCohort," C INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeServiceProvided A on A.PatientSID = C.PatientSID and A.ServiceProvidedCPTSID>0 and A.FeeAuthorizationSID>0 INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeInitialTreatment B on A.FeeAuthorizationSID=B.FeeAuthorizationSID and B.FeeInitialTreatmentSID=A.FeeInitialTreatmentSID and B.Sta3n=A.Sta3n LEFT JOIN MVP_CDWWork.CDWWork.Dim.FeePurposeOfVisit D on A.FeePurposeOfVisitSID=D.FeePurposeOfVisitSID INNER JOIN MVP_CDWWork.CDWWork.Dim.CPT E on A.ServiceProvidedCPTSID=E.CPTSID and E.Sta3n=A.Sta3n WHERE C.C=",cohort),collapse="")
    FeeS=sqlQuery(ch2,query=q,stringsAsFactors=F)
    FeeS=data.table(FeeS)
    FeeS[,InitialTreatmentDateTime:=as.numeric(julian(InitialTreatmentDateTime,origin=as.Date("2000-01-01")))/365.25]
    saveRDS(FeeS,file=paste(c("FeeService.",cohort,".RDS"),collapse=""))
    print(c(dim(FeeS),cohort,date()))
}

for (cohort in seq(0,100)){
    q=paste(c("SELECT C.PatientICN, A.Sta3n, B.InitialTreatmentDateTime, D.ICDCode FROM ",ndbCohort," C INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeServiceProvided A on A.PatientSID = C.PatientSID and A.FeeAuthorizationSID>0 INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeInitialTreatment B on A.FeeAuthorizationSID=B.FeeAuthorizationSID and B.FeeInitialTreatmentSID=A.FeeInitialTreatmentSID and B.Sta3n=A.Sta3n LEFT JOIN MVP_CDWWork.CDWWork.Dim.ICD D on A.ICD9SID=D.ICDSID WHERE A.ICD9SID > 0 and C.C=",cohort),collapse="")
    FeeS=sqlQuery(ch2,query=q,stringsAsFactors=F)
    FeeS=data.table(FeeS)
    FeeS[,InitialTreatmentDateTime:=as.numeric(julian(InitialTreatmentDateTime,origin=as.Date("2000-01-01")))/365.25]
    saveRDS(FeeS,file=paste(c("FeeService.Dx9.",cohort,".RDS"),collapse=""))
    print(c(dim(FeeS),cohort,date()))
}

for (cohort in seq(0,100)){
    q=paste(c("SELECT C.PatientICN, A.Sta3n, B.InitialTreatmentDateTime, D.ICD10Code FROM ",ndbCohort," C INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeServiceProvided A on A.PatientSID = C.PatientSID and A.FeeAuthorizationSID>0 INNER JOIN MVP_CDWWork.CDWWork.Fee.FeeInitialTreatment B on A.FeeAuthorizationSID=B.FeeAuthorizationSID and B.FeeInitialTreatmentSID=A.FeeInitialTreatmentSID and B.Sta3n=A.Sta3n LEFT JOIN MVP_CDWWork.CDWWork.Dim.ICD10 D on A.ICD10SID=D.ICD10SID WHERE A.ICD10SID > 0 and C.C=",cohort),collapse="")
    FeeS=sqlQuery(ch2,query=q,stringsAsFactors=F)
    FeeS=data.table(FeeS)
    FeeS[,InitialTreatmentDateTime:=as.numeric(julian(InitialTreatmentDateTime,origin=as.Date("2000-01-01")))/365.25]
    saveRDS(FeeS,file=paste(c("FeeService.Dx10.",cohort,".RDS"),collapse=""))
    print(c(dim(FeeS),cohort,date()))
}

#Outpatient Procedures
for (cohort in seq(0,100)){ 
    q=paste(c("select C.PatientICN, A.VisitSID, A.sta3n, A.VisitDateTime, A.VProcedureDateTime, B.CPTCode from MVP_CDWWork.CDWWork.Outpat.VProcedure A INNER JOIN MVP_CDWWork.CDWWork.Dim.CPT B on A.CPTSID = B.CPTSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE C.C = ",cohort),collapse="") 
    CPT=sqlQuery(ch2,query=q,stringsAsFactors=F) 
    CPT=data.table(CPT) 
    CPT[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25] 
    CPT[,VProcedureDateTime:=as.numeric(julian(VProcedureDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(CPT,PatientICN) 
    saveRDS(CPT,file=paste(c("oCPT.",cohort,".RDS"),collapse=""))   #this file is in ../tmp, then split below into EM, Med, and other.
    print(c(dim(CPT),cohort,"all"))
} 
for (cohort in seq(0,100)){
    d=readRDS(paste(c("oCPT.",cohort,".RDS"),collapse=""))
    d[,"NumCode":=as.integer(CPTCode)]
    saveRDS(d[NumCode >= 99200 & NumCode < 99500,c(1,2,3,4,5,6)],file=paste(c("oCPT.EM.",cohort,".RDS"),collapse=""))
    saveRDS(d[NumCode >= 90400 & NumCode < 99200,c(1,2,3,4,5,6)],file=paste(c("oCPT.Med.",cohort,".RDS"),collapse=""))
    saveRDS(d[NumCode < 90400 | NumCode > 99500 | is.na(NumCode)==T,c(1,2,3,4,5,6)],file=paste(c("oCPT.Other.",cohort,".RDS"),collapse=""))
    print(c(dim(d),cohort))
}

#outpatient facility attributes
for (cohort in seq(0,100)){
    q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, B.OfficialVAName, B.City, B.Zip, D.StateAbbrev, E.StopCode, E.StopCodeName, F.StopCode, F.StopCodeName, G.FacilityType, G.FullName from MVP_CDWWork.CDWWork.Outpat.Workload A LEFT JOIN MVP_CDWWork.CDWWork.Dim.Institution B on A.InstitutionSID=B.InstitutionSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.State D on B.StateSID=D.StateSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.StopCode E on A.PrimaryStopCodeSID = E.StopCodeSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.StopCode F on A.SecondaryStopCodeSID = F.StopCodeSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.FacilityType G on B.FacilityTypeSID=G.FacilityTypeSID WHERE C.C = ",cohort),collapse="")
    oFac=sqlQuery(ch2,query=q,stringsAsFactors=F) 
    oFac=data.table(oFac) 
    oFac[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(oFac,PatientICN)
    saveRDS(oFac,file=paste(c("oFac.",cohort,".RDS"),collapse="")) 
    print(c(dim(oFac),cohort))
}

#outpatient demographics
for (cohort in seq(0,100)){
    q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, A.ProviderCount, A.DiagnosisCount, A.ProcedureCount, A.PatientFIPS, A.PatientZIP, A.County, A.PatientMaritalStatus, A.PatientReligion, A.PatientNumberOfDependents, A.PatientIncome from MVP_CDWWork.CDWWork.Outpat.Workload A INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE C.C = ",cohort),collapse="")
    oDemo=sqlQuery(ch2,query=q,stringsAsFactors=F) 
    oDemo=data.table(oDemo)
    oDemo[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(oDemo,PatientICN)
    saveRDS(oDemo,file=paste(c("oDemo.",cohort,".RDS"),collapse="")) 
    print(c(dim(oDemo),cohort))
}

#outpatient provider type
for (cohort in seq(0,100)){
    q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, B.ProviderType, B.AreaOfSpecialization from MVP_CDWWork.CDWWork.Outpat.VProvider A LEFT JOIN MVP_CDWWork.CDWWork.Dim.ProviderType B on A.ProviderTypeSID=B.ProviderTypeSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE C.C = ",cohort),collapse="")
    oProv=sqlQuery(ch2,query=q,stringsAsFactors=F) 
    oProv=data.table(oProv)
    oProv[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(oProv,PatientICN)
    saveRDS(oProv,file=paste(c("oProv.",cohort,".RDS"),collapse="")) 
    print(c(dim(oProv),cohort))
}

#outpatient Dx
for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, B.ICDCode, A.PrimarySecondary from MVP_CDWWork.CDWWork.Outpat.VDiagnosis A INNER Join MVP_CDWWork.CDWWork.Dim.ICD B on A.ICD9SID=B.ICDSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE A.ICD9SID != 0 AND C.C = ",cohort),collapse="")
Dx=sqlQuery(ch2,query=q,stringsAsFactors=F)
Dx=data.table(Dx)
Dx[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(Dx,PatientICN)
saveRDS(Dx,file=paste(c("oDx9.",cohort,".RDS"),collapse=""))
print(c(dim(Dx),cohort))
}

for (cohort in seq(0,100)){
q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, B.ICD10Code, A.PrimarySecondary from MVP_CDWWork.CDWWork.Outpat.VDiagnosis A INNER Join MVP_CDWWork.CDWWork.Dim.ICD10 B on A.ICD10SID=B.ICD10SID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE A.ICD10SID != 0 AND C.C = ",cohort),collapse="")
Dx=sqlQuery(ch2,query=q,stringsAsFactors=F)
Dx=data.table(Dx)
Dx[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25]
setkey(Dx,PatientICN)
saveRDS(Dx,file=paste(c("oDx10.",cohort,".RDS"),collapse=""))
print(c(dim(Dx),cohort))
}

#outpatient health factors
for (cohort in seq(0,99)){
    q=paste(c("select C.PatientICN, A.VisitSID, A.VisitDateTime, B.HealthFactorType, A.LevelSeverity, B.HealthFactorCategory from MVP_CDWWork.CDWWork.HF.HealthFactor A INNER JOIN MVP_CDWWork.CDWWork.Dim.HealthFactorType B on A.HealthFactorTypeSID=B.HealthFactorTypeSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE C.C = ",cohort),collapse="")
    HF=sqlQuery(ch2,query=q,stringsAsFactors=F)
    HF=data.table(HF)
    HF[,VisitDateTime:=as.numeric(julian(VisitDateTime,origin=as.Date("2000-01-01")))/365.25]
    setkey(HF,PatientICN)
    saveRDS(HF,file=paste(c("HF.",cohort,".RDS"),collapse=""))
    print(c(dim(HF),cohort))
}

q='select BenMcMahonTemp.kdi\mcmahonb.mvp011w.cc022'
#major CPT categories (not used now when pulling data)
#",ndbCohort," C
#         <NA>
#2                           *Missing*
#3              *Unknown at this time*
#4                          ANESTHESIA
#5           EVALUATION AND MANAGEMENT
#6  EVALUATION AND MANAGEMENT SERVICES
#7                          LABORATORY
#8                            MEDICINE
#9            PATHOLOGY AND LABORATORY
#10                          RADIOLOGY
#11                            SURGERY

#laboratory test results.  For the first 34 cohorts, data are pulled and merged in one table.  for cohorts 35-100, each cohort is subset into 100 subgroups for both data pulling and merging / splitting, both for memory reasons and to enable parallel database queries.

for (cohort in seq(78,60.-1)){ q=paste(c("select C.PatientICN, A.LabChemTestSID, A.LabChemSpecimenDateTime, B.LabChemTestName, A.LabPanelSID, E.Topography, B.LabTestType, D.LOINC, A.Sta3n, A.LabChemResultValue, A.LabChemResultNumericValue, A.Units, A.Abnormal, A.RefHigh, A.RefLow from MVP_CDWWork.CDWWork.Chem.PatientLabChem A INNER JOIN MVP_CDWWork.CDWWork.Dim.LabChemTest B on A.LabChemTestSID = B.LabChemTestSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.LOINC D on A.LOINCSID = D.LOINCSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.Topography E on A.TopographySID=E.TopographySID WHERE C.C = ",cohort),collapse="") 
labs=sqlQuery(ch2,query=q,stringsAsFactors=F) 
labs=data.table(labs) 
labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25] 
setkey(labs,PatientICN) 
saveRDS(labs,file=paste(c("labs.",cohort,".RDS"),collapse="")) 
print(c(dim(labs),cohort)) } 

LL=readRDS("Loinc.labs.RDS")  # Ben's top 4,000 LOINC codes, classified hierarchically to split data
LL3=LL[,c(1,3,15,16)]
LL3$class=as.character(LL3$class)
LL3$subclass=as.character(LL3$subclass)
colnames(LL3)[3]="Lclass"
colnames(LL3)[4]="Lsubclass"
setkey(LL3,LOINC)
MVPLabs=readRDS("DimLab18.RDS")          # Kelly's team's classification of 110,719 LabChemTestSIDs, (25,485 LabChemTestNames) mapped to 108 short lab names
MVP=unique(MVPLabs,by=c("LabChemTestSID","Sta3n","Topography"))
MVQ=MVP[,c(2,4,7,10,9)]
SNmap=readRDS("MVPLabs.ShortName.Hier.RDS")
setkey(SNmap,ShortName)
setkey(MVPLabs,ShortName)
x=SNmap[MVPLabs]
MVP=unique(x[,c(6,3,4,1)],by="LabChemTestSID")
setkey(MVP,LabChemTestSID)

for (cohort in seq(0,30)){
    labs=readRDS(file=paste(c("labs.",cohort,".RDS"),collapse=""))
    setkey(labs,LabChemTestSID)    
    labs=MVP[labs]
    setkey(labs,LOINC)
    setkey(LL3,LOINC)
    labs=LL3[labs]
    labs[which(is.na(class))]$class=labs[which(is.na(class)),]$Lclass
    labs[which(is.na(subclass))]$subclass=labs[which(is.na(subclass)),]$Lsubclass
    labs[,Lclass:=NULL]
    labs[,Lsubclass:=NULL]
    labs=labs[,c(7,8,3,9,6,1,4,5,10,11,12,13,14,15,16,17,18,19)]
    setkey(labs,PatientICN)
        labs=MVQ[labs,on=c("LabChemTestSID","Sta3n","Topography")]
        if (!is.numeric(labs[1,LabChemSpecimenDateTime])) labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25]
        labs[,i.ShortName:=NULL]
        labs=labs[,c(6,7,2,1,8,3,4,5,9,10,11,12,13,14,15,16,17,18,19)]
        setkey(labs,PatientICN)
     #   saveRDS(labs,file=paste(c("labout/",ft,cohort,".RDS"),collapse=""))
    saveRDS(labs[which(is.na(class)),],file=paste(c("labout/Other.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="M" & subclass=="electrolytes",],file=paste(c("labout/M.electrolytes.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="M" & (subclass=="kidney" | subclass=="protein"),],file=paste(c("labout/M.kidney.protein.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="M" & !(subclass=="kidney"|subclass=="protein"|subclass=="electrolytes"),],file=paste(c("labout/M.protein.urine.vital.o.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="H" & (subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"),],file=paste(c("labout/H.WBC.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="H" & (subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("labout/H.RBC.iron.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="H" & !(subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"|subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("labout/H.CBC.coag.o.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="L" | class=="S",],file=paste(c("labout/Lipid.Sugars.",cohort,".RDS"),collapse=""))
    saveRDS(labs[class=="E" | class=="D",],file=paste(c("labout/Enzymes.hormones.Drugs.",cohort,".RDS"),collapse=""))
    saveRDS(labs[!(class=="M"|class=="H"|class=="H"|class=="L"|class=="S"|class=="E"|class=="D"),],file=paste(c("labout/Pathogens.Cancer.GIAx.",cohort,".RDS"),collapse=""))
}

CCa=readRDS("Cohort.2022.RDS")
#CC=unique(C,by="PatientSID")
#CC$PatientSID=as.integer("PatientSID")
#CCa=CC[!is.na(C),]

setkey(CCa,PatientSID)

for (cohort in seq(49,49)){   #for the DoB cohorts with larger number of labs, subset each by 100 so data pulls can be done in parallel without running out of memory in R session, typically ~5 jobs per VM.  database server has 20 cores devoted to pulling data.
    CCb=CCa[C==cohort,] 
    setkey(CCb,PatientSID) 
    N=dim(CCb)[1]  
    for (ss in seq(27,99)){ 
        #mnYoB=min(CC[C==cohort,YoB]) #mxYoB=max(CC[C==cohort,YoB]) #mYoB=median(CC[C==cohort,YoB]) 
        a=CCb$PatientSID[ss*N/100+1] 
        b=CCb$PatientSID[(ss+1)*N/100] 
        q=paste(c("select C.PatientICN, A.LabChemTestSID, A.LabChemSpecimenDateTime, B.LabChemTestName, A.LabPanelSID, E.Topography, B.LabTestType, D.LOINC, A.Sta3n, A.LabChemResultValue, A.LabChemResultNumericValue, A.Units, A.Abnormal, A.RefHigh, A.RefLow from MVP_CDWWork.CDWWork.Chem.PatientLabChem A INNER JOIN MVP_CDWWork.CDWWork.Dim.LabChemTest B on A.LabChemTestSID = B.LabChemTestSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.LOINC D on A.LOINCSID = D.LOINCSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.Topography E on A.TopographySID=E.TopographySID WHERE A.PatientSID BETWEEN ",a," AND ",b," AND C.C = ",cohort),collapse="") 
        date() 
        labs=sqlQuery(ch2,query=q,stringsAsFactors=F) 
        labs=data.table(labs) 
        labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25] 
        setkey(labs,PatientICN) 
        saveRDS(labs,file=paste(c("labs/labs.",cohort,"a.",ss,".RDS"),collapse="")) 
        print(c(dim(labs),cohort,"a",ss,date())) 
    }
}

LL=readRDS("Loinc.labs.RDS")  # Ben's top 4,000 LOINC codes, classified hierarchically to split data
LL3=LL[,c(1,3,15,16)]
LL3$class=as.character(LL3$class)
LL3$subclass=as.character(LL3$subclass)
colnames(LL3)[3]="Lclass"
colnames(LL3)[4]="Lsubclass"
setkey(LL3,LOINC)
MVPLabs=readRDS("DimLab18.RDS")          # Kelly's team's classification of 110,719 LabChemTestSIDs, (25,485 LabChemTestNames) mapped to 108 short lab names
MVP=unique(MVPLabs,by=c("LabChemTestSID","Sta3n","Topography"))
MVQ=MVP[,c(2,4,7,10,9)]
SNmap=readRDS("MVPLabs.ShortName.Hier.RDS")
setkey(SNmap,ShortName)
setkey(MVPLabs,ShortName)
x=SNmap[MVPLabs]
MVP=unique(x[,c(6,3,4,1)],by="LabChemTestSID")
setkey(MVP,LabChemTestSID)

for (cohort in seq(34,100,)){
labs=1;rm(labs)
  flag=0
 for (i in seq(0,49)){
    labsi=readRDS(file=paste(c("labs.tmp/labs.",cohort,"a.",i,".RDS"),collapse=""))
    setkey(labsi,LabChemTestSID)    
    labsi=MVP[labsi]
    setkey(labsi,LOINC)
    setkey(LL3,LOINC)
    labsi=LL3[labsi]
    labsi[which(is.na(class))]$class=labsi[which(is.na(class)),]$Lclass
    labsi[which(is.na(subclass))]$subclass=labsi[which(is.na(subclass)),]$Lsubclass
    labsi[,Lclass:=NULL]
    labsi[,Lsubclass:=NULL]
    labsi=labsi[,c(7,8,3,9,6,1,4,5,10,11,12,13,14,15,16,17,18,19)]
    if (flag==0){labs=labsi} else {labs=rbind(labs,labsi)}
  flag=1
  cat("read",cohort,".",i," size is ",dim(labsi),"\n")
 }
        labs=MVQ[labs,on=c("LabChemTestSID","Sta3n","Topography")]
        if (!is.numeric(labs[1,LabChemSpecimenDateTime])) labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25]
        labs[,i.ShortName:=NULL]
        labs=labs[,c(6,7,2,1,8,3,4,5,9,10,11,12,13,14,15,16,17,18,19)]
        setkey(labs,PatientICN)
     #   saveRDS(labs,file=paste(c("labout/",ft,cohort,".RDS"),collapse=""))
 setkey(labs,PatientICN)
#MAKE sure to send missing into MissingLOINC.cohort.RDS!!!
saveRDS(labs[which(is.na(class)),],file=paste(c("Other.",cohort,"a.RDS"),collapse=""))
#labs=LL3[labs]
saveRDS(labs[class=="M" & subclass=="electrolytes",],file=paste(c("M.electrolytes.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="M" & (subclass=="kidney" | subclass=="protein"),],file=paste(c("M.kidney.protein.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="M" & !(subclass=="kidney"|subclass=="protein"|subclass=="electrolytes"),],file=paste(c("M.protein.urine.vital.o.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="H" & (subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"),],file=paste(c("H.WBC.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="H" & (subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("H.RBC.iron.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="H" & !(subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"|subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("H.CBC.coag.o.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="L" | class=="S",],file=paste(c("Lipid.Sugars.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[class=="E" | class=="D",],file=paste(c("Enzymes.hormones.Drugs.",cohort,"a.RDS"),collapse=""))
saveRDS(labs[!(class=="M"|class=="H"|class=="H"|class=="L"|class=="S"|class=="E"|class=="D"),],file=paste(c("Pathogens.Cancer.GIAx.",cohort,"a.RDS"),collapse=""))

labs=1;rm(labs)
  flag=0
 for (i in seq(50,99)){
    labsi=readRDS(file=paste(c("labs.tmp/labs.",cohort,"a.",i,".RDS"),collapse=""))
    setkey(labsi,LabChemTestSID)    
    labsi=MVP[labsi]
    setkey(labsi,LOINC)
    setkey(LL3,LOINC)
    labsi=LL3[labsi]
    labsi[which(is.na(class))]$class=labsi[which(is.na(class)),]$Lclass
    labsi[which(is.na(subclass))]$subclass=labsi[which(is.na(subclass)),]$Lsubclass
    labsi[,Lclass:=NULL]
    labsi[,Lsubclass:=NULL]
    labsi=labsi[,c(7,8,3,9,6,1,4,5,10,11,12,13,14,15,16,17,18,19)]
    if (flag==0){labs=labsi} else {labs=rbind(labs,labsi)}
  flag=1
  cat("read",cohort,".",i," size is ",dim(labsi),"\n")
 }
        labs=MVQ[labs,on=c("LabChemTestSID","Sta3n","Topography")]
        if (!is.numeric(labs[1,LabChemSpecimenDateTime])) labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25]
        labs[,i.ShortName:=NULL]
        labs=labs[,c(6,7,2,1,8,3,4,5,9,10,11,12,13,14,15,16,17,18,19)]
        setkey(labs,PatientICN)
     #   saveRDS(labs,file=paste(c("labout/",ft,cohort,".RDS"),collapse=""))
 setkey(labs,PatientICN)
#MAKE sure to send missing into MissingLOINC.cohort.RDS!!!
saveRDS(labs[which(is.na(class)),],file=paste(c("Other.",cohort,"b.RDS"),collapse=""))
#labs=LL3[labs]
saveRDS(labs[class=="M" & subclass=="electrolytes",],file=paste(c("M.electrolytes.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="M" & (subclass=="kidney" | subclass=="protein"),],file=paste(c("M.kidney.protein.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="M" & !(subclass=="kidney"|subclass=="protein"|subclass=="electrolytes"),],file=paste(c("M.protein.urine.vital.o.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="H" & (subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"),],file=paste(c("H.WBC.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="H" & (subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("H.RBC.iron.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="H" & !(subclass=="WBC"|subclass=="diff WBC"|subclass=="flow"|subclass=="RBC"|subclass=="hct/Hb"|subclass=="heme"|subclass=="iron"),],file=paste(c("H.CBC.coag.o.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="L" | class=="S",],file=paste(c("Lipid.Sugars.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[class=="E" | class=="D",],file=paste(c("Enzymes.hormones.Drugs.",cohort,"b.RDS"),collapse=""))
saveRDS(labs[!(class=="M"|class=="H"|class=="H"|class=="L"|class=="S"|class=="E"|class=="D"),],file=paste(c("Pathogens.Cancer.GIAx.",cohort,"b.RDS"),collapse=""))
}

#merge files a and b, split above for memory reasons (64GB/VM)
for (cohort in seq(50,100)){
    for (ft in c("Other.","M.electrolytes.","M.kidney.protein.","M.protein.urine.vital.o.","H.WBC.","H.RBC.iron.","H.CBC.coag.o.","Lipid.Sugars.","Enzymes.hormones.Drugs.","Pathogens.Cancer.GIAx.")){
        labs=readRDS(paste(c(ft,cohort,"a.RDS"),collapse=""))
        labsb=readRDS(paste(c(ft,cohort,"b.RDS"),collapse=""))
        labs=rbind(labs,labsb)
        setkey(labs,PatientICN)
        saveRDS(labs,file=paste(c("labout/",ft,cohort,".RDS"),collapse=""))
        cat("wrote",cohort,ft,"\n")
    }
}

#re-join MVPShortnames with LabChemTestSID + Sta3n + Topology; change dates to nuneric, if needed
MVPLabs=readRDS("DimLab18.RDS")          # Kelly's team's classification of 110,719 LabChemTestSIDs, (25,485 LabChemTestNames) mapped to 108 short lab names
MVP=unique(MVPLabs,by=c("LabChemTestSID","Sta3n","Topography"))
MVQ=MVP[,c(2,4,7,10,9)]
for (cohort in seq(0,100)){
    for (ft in c("Other.","M.electrolytes.","M.kidney.protein.","M.protein.urine.vital.o.","H.WBC.","H.RBC.iron.","H.CBC.coag.o.","Lipid.Sugars.","Enzymes.hormones.Drugs.","Pathogens.Cancer.GIAx.")){
        d=readRDS(paste(c(ft,cohort,".RDS"),collapse=""))
        e=MVQ[d,on=c("LabChemTestSID","Sta3n","Topography")]
        if (!is.numeric(e[1,LabChemSpecimenDateTime])) e[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25]
        e[,i.ShortName:=NULL]
        e=e[,c(6,7,2,1,8,3,4,5,9,10,11,12,13,14,15,16,17,18,19)]
        setkey(e,PatientICN)
        saveRDS(e,file=paste(c("labout/",ft,cohort,".RDS"),collapse=""))
    }
}
#
#breakdown:of Ben's classification, on MVP Cohort 17 or 18.
 #            M (metabolites)      58,024,428  = M.electrolytes, 22M  + M.kidney/stone/protein, 19M  + M.blank,neurotrans, urine,note,stool,vital, 17M

#             H (hematology, coag) 58,997,619  = H.WBC,diffWBC,flow, 26M  + H.RBC,iron,heme,hct/Hb, 20M + H.blank,CBC,coag,calcitonin,  13M

#             L (lipids)            7,489,648
#             S (sugars)            9,411,590

#             E (enzyme, hormone)  10,299,186
#             D (drug, vitamin)     5,321,264

#             C (cancer)            1,370,313
#             G (genetic, reprod)      62,436
#             P (pathogens)         1,722,293
#             I (autoimmune)          228,833
#             A (allergies)            60,656
#             X (?)                     1,850
#             Missing:                  1,574


#splitting in half (using mYoB, above), and getting those specifically with missing loinc, although that should already be included, above
#for (cohort in c(78,82,84,86,88)){
#q=paste(c("select C.PatientICN, A.LabChemTestSID, A.LabChemSpecimenDateTime, B.LabChemTestName, B.LabChemPrintTestName, B.LabTestType, D.LOINC, A.Sta3n, A.LabChemResultValue, A.LabChemResultNumericValue, A.Units, A.Abnormal, A.RefHigh, A.RefLow from MVP_CDWWork.CDWWork.Chem.PatientLabChem A INNER JOIN MVP_CDWWork.CDWWork.Dim.LabChemTest B on A.LabChemTestSID = B.LabChemTestSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID INNER JOIN MVP_CDWWork.CDWWork.Dim.LOINC D on A.LOINCSID = D.LOINCSID WHERE D.LOINC = '*Missing*' AND C.C = ",cohort),collapse="")
#date()
#labs=sqlQuery(ch2,query=q,stringsAsFactors=F)
#labs=data.table(labs)
#labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25] 
#setkey(labs,PatientICN)
#saveRDS(labs,file=paste(c("MissingLOINC.",cohort,".RDS"),collapse=""))
#print(c(dim(labs),cohort," a ",date()))
#}

#q=paste(c("select C.PatientICN, A.LabChemTestSID, A.LabChemSpecimenDateTime, B.LabChemTestName, B.LabChemPrintTestName, B.LabTestType, D.LOINC, A.Sta3n, A.LabChemResultValue, A.LabChemResultNumericValue, A.Units, A.Abnormal, A.RefHigh, A.RefLow from MVP_CDWWork.CDWWork.Chem.PatientLabChem A INNER JOIN MVP_CDWWork.CDWWork.Dim.LabChemTest B on A.LabChemTestSID = B.LabChemTestSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID INNER JOIN MVP_CDWWork.CDWWork.Dim.LOINC D on A.LOINCSID = D.LOINCSID WHERE C.YoB >= ",mYoB," AND C.C = ",cohort),collapse="")
#labs=sqlQuery(ch2,query=q,stringsAsFactors=F)
#labs=data.table(labs)
#labs[,LabChemSpecimenDateTime:=as.numeric(julian(LabChemSpecimenDateTime,origin=as.Date("2000-01-01")))/365.25] 
#setkey(labs,PatientICN)
#saveRDS(labs,file=paste(c("labs.",cohort,"b.RDS"),collapse=""))
#print(c(dim(labs),cohort,"b"))
#}

#MH surveys
for (cohort in seq(0,99)){
    q=paste(c("select C.PatientICN, A.SurveyName, A.SurveyGivenDateTime, A.SurveyScale, A.RawScore, A.TransformedScore1 from MVP_CDWWork.CDWWork.MH.SurveyResult A INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID WHERE C.C = ",cohort),collapse="")
    MH=sqlQuery(ch2,query=q,stringsAsFactors=F)
    MH=data.table(MH)
    MH[,SurveyGivenDateTime:=as.numeric(julian(SurveyGivenDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(MH,PatientICN)
    saveRDS(MH,file=paste(c("SResult.",cohort,".RDS"),collapse=""))
    print(c(dim(MH),cohort))
}

for (cohort in seq(0,99)){
    q=paste(c("select C.PatientICN, A.SurveyName, A.SurveyGivenDateTime, B.SurveyQuestionText, D.SurveyChoiceText from MVP_CDWWork.CDWWork.MH.SurveyAnswer A LEFT JOIN MVP_CDWWork.CDWWork.Dim.SurveyQuestion B on A.SurveyQuestionSID=B.SurveyQuestionSID INNER JOIN ",ndbCohort," C on A.PatientSID = C.PatientSID LEFT JOIN MVP_CDWWork.CDWWork.Dim.SurveyChoice D on A.SurveyChoiceSID=D.SurveyChoiceSID WHERE C.C = ",cohort),collapse="")
    MH=sqlQuery(ch2,query=q,stringsAsFactors=F)
    MH=data.table(MH)
    MH[,SurveyGivenDateTime:=as.numeric(julian(SurveyGivenDateTime,origin=as.Date("2000-01-01")))/365.25] 
    setkey(MH,PatientICN)
    saveRDS(MH,file=paste(c("SQuestions.",cohort,".RDS"),collapse=""))
    print(c(dim(MH),cohort))
}

#compute median bp, ht, wt, BMI, date of first and most recent (vital or pharmacy or CPT), number of Years, Quarters, months, days with (vital or pharmacy or CPT)
for (cohort in seq(0,100)){
d=readRDS(paste(c("../tmp/DemoRE.",cohort,".RDS"),collapse=""))
bp=readRDS(paste(c("vital/Vital.bp.",cohort,".RDS"),collapse=""))
Vo=readRDS(paste(c("vital/Vital.other.",cohort,".RDS"),collapse=""))
CN=readRDS(paste(c("meds/Rx.CN.",cohort,".RDS"),collapse=""))
CV=readRDS(paste(c("meds/Rx.CN.",cohort,".RDS"),collapse=""))
A=readRDS(paste(c("meds/Rx.CN.",cohort,".RDS"),collapse=""))
HS=readRDS(paste(c("meds/Rx.CN.",cohort,".RDS"),collapse=""))
Oth=readRDS(paste(c("meds/Rx.CN.",cohort,".RDS"),collapse=""))
EM=readRDS(paste(c("outpat/oCPT.EM.",cohort,".RDS"),collapse=""))
O=readRDS(paste(c("outpat/oCPT.Other.",cohort,".RDS"),collapse=""))
Med=readRDS(paste(c("outpat/oCPT.Med.",cohort,".RDS"),collapse=""))

bpM=bp[,median(Systolic,by="PatientICN")]

uRx=rbind(CN[,c("PatientICN","FillDateTime")],CV[,c("PatientICN","FillDateTime")],HS[,c("PatientICN","FillDateTime")],A[,c("PatientICN","FillDateTime")],Oth[,c("PatientICN","FillDateTime")])
uV=rbind(bp[,c("PatientICN","VitalSignTakenDateTime")],Vo[,c("PatientICN","VitalSignTakenDateTime")])
uCPT=rbind(EM[,c("PatientICN","VisitDateTime")],Med[,c("PatientICN","VisitDateTime")],O[,c("PatientICN","VisitDateTime")])
colnames(uRx)[2]="Date"
colnames(uV)[2]="Date"
colnames(uCPT)[2]="Date"
u=rbind(uRx,uV,uCPT)
uF=u[,min(Date),by="PatientICN"]
uL=u[,max(Date),by="PatientICN"]
colnames(uF)[2]="FirstDate"
colnames(uL)[2]="LastDate"
#u[,"S":=substr(as.Date(Date*365.25,origin=as.Date("2000-01-01")),1,4)]
u[,"S":=floor(Date)]
bins=unique(u,by=c("PatientICN","S"))
uY=bins[Date>0,.N,by="PatientICN"]
colnames(uY)[2]="NYears"
u[,"S":=floor(Date*4)]
bins=unique(u,by=c("PatientICN","S"))
uQ=bins[Date>0,.N,by="PatientICN"]
colnames(uQ)[2]="NQuarters"
#u[,"S":=substr(as.Date(Date*365.25,origin=as.Date("2000-01-01")),1,7)]
u[,"S":=floor(Date*12)]
bins=unique(u,by=c("PatientICN","S"))
uM=bins[Date>0,.N,by="PatientICN"]
colnames(uM)[2]="NMonths"
#u[,"S":=substr(as.Date(Date*365.25,origin=as.Date("2000-01-01")),1,10)]
u[,"S":=floor(Date*365)]
bins=unique(u,by=c("PatientICN","S"))
uD=bins[Date>0,.N,by="PatientICN"]
colnames(uD)[2]="NDays"

Ht=Vo[VitalTypeAbbreviation=="HT" & VitalResultNumeric>36,median(VitalResultNumeric),by="PatientICN"]
colnames(Ht)[2]="Height"
Wt=Vo[VitalTypeAbbreviation=="WT" & VitalResultNumeric>50,median(VitalResultNumeric),by="PatientICN"]
colnames(Wt)[2]="Weight"
Systolic=bp[,median(Systolic),by="PatientICN"]
colnames(Systolic)[2]="Systolic"

util=list(d,Ht,Wt,Systolic,uF,uL,uY,uQ,uM,uD)
demoOut=Reduce(function(...) merge(...,all=T,by="PatientICN"),util)
demoOut[,"BMI":=703*Weight/Height/Height]
demo=demoOut[,c(seq(1,10),18,seq(11,17))]
setkey(demo,PatientICN)
saveRDS(demo,file=paste(c("DemoU.",cohort,".RDS"),collapse=""))
print(c(dim(demo),cohort,date()))
}