import os
if os.getcwd().endswith('Setup'):
    os.chdir('..')

import NotebookImport
from Setup.Imports import *

"""
Explanation: Read HIV Data
End of explanation
"""


c1 = pd.read_excel(ucsd_path + 'DESIGN_Fox_v2_Samples-ChipLAyout-Clinical UNMC-UCSD methylomestudy.xlsx', 
                   'HIV- samples from OldStudy', index_col=0)
c2 = pd.read_excel(ucsd_path + 'DESIGN_Fox_v2_Samples-ChipLAyout-Clinical UNMC-UCSD methylomestudy.xlsx', 
                   'HIV+ samples', index_col=0)
clinical = c1.append(c2)
clinical['Sentrix_Position'] = clinical['Sentrix_Position\\'].map(lambda s: s[:-1])
del clinical['Sentrix_Position\\']

"""
Explanation: Read in Clinical Data
End of explanation
"""


age_new = pd.read_csv(ucsd_path + 'UpdatesAges-Infection.csv', index_col=0)
age = age_new.age.combine_first(clinical.age)
age.name= 'age'
clinical['age'] = age
l = 'estimated duration hiv (months)'
clinical[l] = age_new['Estimated Duration HIV+ (months)'].combine_first(clinical[l])

"""
Explanation: Update clinical data with new data provided by Howard Fox
End of explanation
"""


diabetes = clinical['diabetes'].combine_first(clinical['Diabetes @ 000'])
diabetes = diabetes.replace('N','no')
clinical['diabetes'] = diabetes
del clinical['Diabetes @ 000']
diabetes.value_counts()

"""
Explanation: Clean up diabetes across annotation files
End of explanation
"""


ethnicity = clinical.ethnicity
ethnicity = ethnicity.replace('wht','white')
ethnicity = ethnicity.replace('Caucasian - European','white')
clinical['ethnicity'] = ethnicity
ethnicity.value_counts()

"""
Explanation: All of the patients are white or Caucasian
End of explanation
"""


clinical['sex'] = clinical['sex'].fillna('M')

"""
Explanation: Sex is not recorded for the cases but they are all HIV+ men.
End of explanation
"""


bmi = clinical['bmi'].combine_first(clinical['BMI'])
clinical['BMI'] = bmi
clinical = clinical[clinical.columns.diff(['bmi'])]

#Do not import
bmi.hist()

current_usage = ["Current 'Other' dx", 'Current Alcohol dx',
                 'Current Bipolar I', 'Current Bipolar II',
                 'Current Cannabis dx', 'Current Cocaine dx',
                 'Current Dysthymia', 'Current Halucinogen dx',
                 'Current Inhalant dx', 'Current MDD',
                 'Current Methamphetamine dx', 'Current Opioid dx',
                 'Current PCP dx', 'Current Sedative dx',
                 'Any Current Substance dx']
current_usage = clinical[current_usage]
current_usage.dropna(how='all').apply(pd.value_counts).fillna(0).T

past_usage = ["LT 'Other' dx",  'LT Alcohol dx', 'LT Bipolar I', 
              'LT Bipolar II', 'LT Cannabis dx', 'LT Cocaine dx',
              'LT Dysthymia', 'LT Halucinogen dx', 'LT Inhalant dx',
              'LT MDD', 'LT Methamphetamine dx', 'LT Opioid dx',
              'LT PCP dx', 'LT Sedative dx', 'Any LT Substance dx']
past_usage = clinical[past_usage]
past_usage.dropna(how='all').apply(pd.value_counts).fillna(0).T

usage = current_usage.join(past_usage).dropna(how='all')
usage.to_csv(FIGDIR + 'drug_usage.csv')

merge_col = ['age','BMI','sex','diabetes','ethnicity','Sample_Plate','Sample_Well',
             'Sentrix_ID','Sentrix_Position']
shared_clinical = clinical[merge_col]
shared_clinical.to_csv(FIGDIR + 'shared_clinical.csv')

c = [u'Site', u'ARV History', u'ARV Status', u'BDI > 17', u'CDC stage', 
     u'estimated duration hiv (months)', u'Current Regimen',
     u'Regimen Type', u'adherence %',
     u'HCV', u'IADL', u'RPR', u'Utox', 
     u'beck total', u'global impairment', u'paofi total',]
hiv_clinical = clinical[c].dropna(axis=0)
hiv_clinical.to_csv(FIGDIR + 'hiv_clinical.csv')

"""
Explanation: Fix BMI to unified labels
End of explanation
"""


clinical.diabetes.value_counts()

clinical['Diabetes @ 000'].value_counts()

"""
Explanation: Trimming the clinical dataset
None of the patients are diabetic
End of explanation
"""


clinical['HCV'].dropna(0).value_counts(0)

"""
Explanation: All of the patients are hepatitis C negative
End of explanation
"""


clinical['ARV History'].value_counts()

clinical['ARV Status'].value_counts()

"""
Explanation: All patients are currently using anti-retoviral therepy, but 5 patients treatment is not classified as HAART
End of explanation
"""


clinical['adherent'].value_counts()

#Do not import
clinical['adherence %'].hist()

"""
Explanation: All patients are reported as adhererent, but a few are not 100% adherent
End of explanation
"""


reg = clinical['Current Regimen'].dropna().str.split('/').map(sorted)
drugs = {r for s in reg for r in s}
drug_mat = pd.DataFrame({i: {d: d in s for d in drugs} for i,s in 
                         reg.iteritems()}).T
drug_mat.sum().order()

"""
Explanation: We have a wide varienty of regimens with 81 unique combinations of 35 drugs
End of explanation
"""


clinical['Regimen Type'].value_counts()

kill_list = ['zhang id', 'diabetes', 'Methylation ID', 'Diabetes @ 000',
             'Sentrix_ID','Sample_Plate','Sample_Well','Sentrix_Position\\',
             ]
drugs = ['ARV History', 'ARV Status', 'Current Regimen', 'Regimen Type',
         'adherence %' ,'adherent']

left = [c for c in clinical if c not in past_usage and c not in current_usage
        and c not in drugs]

age = clinical.age

clinical['BDI > 17'].value_counts()

clinical['ARV History'].value_counts()

clinical['CDC stage'].value_counts()

"""
Explanation: These can be broken down into 8 regimen types
End of explanation
"""


iadl = clinical.IADL
iadl.value_counts()

"""
Explanation: IADL = Instrumental activities of daily living
End of explanation
"""


clinical['global impairment'].value_counts()

"""
Explanation: Global imparement
End of explanation
"""


#Do not import
paofi = clinical['paofi total']
paofi.hist()

"""
Explanation: Patient's Assessment of Own Functioning Inventory (PAOFI)
End of explanation
"""


clinical.RPR.value_counts()

"""
Explanation: RPR (Rapid plasma reagin) is a diagnostic used to detect syphilis
End of explanation
"""


site = clinical.Site
site.value_counts()

clinical.Utox.value_counts()

"""
Explanation: Batches
End of explanation
"""


#Do not import
beck = clinical['beck total'].dropna()
beck.hist()

#Do not import
clinical['paofi total'].hist()

"""
Explanation: Beck Depression Inventory is a questionarre measuring depression levels (from Wikipedia)
* 0–13: minimal depression 
* 14–19: mild depression
 20–28: moderate depression
 29–63: severe depression.
End of explanation
"""


labs = pd.read_excel(ucsd_path + 'fox_methylation_labdata_073014.xlsx',
                     index_col=0)
labs['nb exam year']= labs['nb exam date'].map(lambda s: s.year)
del labs['nb exam date']
labs = labs.dropna(axis=1, how='all')
labs = labs.ix[labs.index.intersection(clinical.index)]

hiv_rna = labs[['rnvalue PLASMA','LLQ PLASMA']]
hiv_rna.to_csv(FIGDIR + 'hiv_rna.csv')

labs.WBC.hist()

chem = ['WBC','RBC','HGB','HCT','MCV',
        'MCH','MCHC','Platelets']
labs[chem].to_csv(FIGDIR + 'blood_workup.csv')

imm= ['CD4 Absolute',
      'CD8 Absolute',
      'CD3 Absolute',
      'CD4/CD8 ratio',
      'Neutrophil %','Lymphocyte %',
      'Monocyte %','Eosinophil %','Basophil %']
cell_types = ['Neutrophil %','Lymphocyte %','Monocyte %',
              'Eosinophil %','Basophil %']
hiv_imm = labs[imm].dropna(how='all')

"""
Explanation: Read in lab blood work data

abstract exam date onto exam year for anonymity
End of explanation
"""


keepers = labs.index.diff(['RG065','RG175','RG279','RA182','RM285'])
hiv_imm[cell_types] = hiv_imm[cell_types].div(hiv_imm[cell_types].sum(1), axis=0) * 100
keep = pd.Series('Keep', keepers).ix[hiv_imm.index].fillna('Reject')
hiv_imm['QC status'] = keep
hiv_imm.to_csv(FIGDIR + 'HIV_cell_comp.csv')

#Do not import
fig, axs = subplots(1,2, figsize=(9,4))
clinical.age.hist(ax=axs[0])
clinical['estimated duration hiv (months)'].hist(ax=axs[1])
axs[0].set_xlabel('Age')
axs[1].set_xlabel('estimated duration hiv (months)')
for ax in axs:
    ax.set_ylabel('# of Patients')
    prettify_ax(ax)

"""
Explanation: Dropping five patients because they don't look Kosher  
Renormalizing cell percentages because some don't sum to 100%
End of explanation
"""


duration = clinical['estimated duration hiv (months)']
duration = (1.*duration.notnull()) + (1.*duration > 100)
duration = duration.map({0:'Control',1:'HIV Short',2:'HIV Long'})
duration.value_counts()

"""
Explanation: As can be seen above, we have two groups of patients with respect to HIV duration, we don't really have the sample size to tease apart any differences other than this main distiction so for now I am just treating duration of HIV as a categorical variable (e.g. controls, short exposure and long exposure)
End of explanation
"""


duration = duration.ix[ti(clinical['ARV Status'] != 'non-HAART')]
duration = duration.ix[ti(clinical.sex != 'F')].dropna()
duration.value_counts()

"""
Explanation: Patient Selection Criteria

There are a couple of female patients in the controls, we are going to get rid of those as all of the cases are males  
Most of the HIV patients are under HAART therepy, there are a few that are not and we are going to filter those out for now and possibly look at them after the primary analysis
End of explanation
"""


df_hiv = pd.read_hdf('/data/methylation_norm.h5', 'quant_BMIQ_adj')
df_hiv = df_hiv.ix[:, duration.index]
df_hiv = df_hiv.dropna(1)

"""
Explanation: Read in HIV Methylation data
Read in quantile-normalized data, adjusted for cellular compositions and then normalized agin using BMIQ.
End of explanation
"""


df_hiv_n = pd.read_hdf('/data/methylation_norm.h5', 'BMIQ_Horvath')
df_hiv_n = df_hiv_n.ix[:, duration.index]
df_hiv_n = df_hiv_n.dropna(1)

"""
Explanation: Read in data processed with BMIQ using Horvath's gold standard.
End of explanation
"""


flow_sorted_data = pd.read_hdf('/data/methylation_annotation.h5','flow_sorted_data_horvath_norm')
cell_type = pd.read_hdf('/data/methylation_annotation.h5', 'label_map')
cell_counts = pd.read_hdf('/data/dx_methylation.h5', 'cell_counts')

n2 = flow_sorted_data.groupby(cell_type, axis=1).mean()
avg = n2[cell_counts.columns].dot(cell_counts.ix[df_hiv.columns].T)
d2 = df_hiv_n.ix[avg.index, df_hiv.columns].dropna(axis=[0,1], how='all')
cc = avg.ix[:, ti(duration=='Control')].mean(1)
df_hiv_n = (d2 - avg).add(cc, axis=0).dropna(how='all')

keepers = duration.index.intersection(df_hiv.columns.intersection(df_hiv_n.columns))
duration = duration.ix[keepers]

duration.value_counts()

store = pd.HDFStore('/data/dx_methylation.h5')
study = store['study']
age = store['age']
gender = store['gender']

"""
Explanation: Adjust this data for cellular composition.  This is done after the normalization to not mess around with Horvath's pipeline too much.
End of explanation
"""


detection_p = pd.read_hdf('/data/dx_methylation.h5', 'detection_p')
#detection_p = detection_p[detection_p[0] > 10e-5]
detection_p = detection_p[detection_p.Sample_Name.isin(duration.index)]

ff = detection_p.groupby('level_0').size() > 3
ff.value_counts()

STORE = '/data/methylation_annotation.h5'
snps = pd.read_hdf(STORE, 'snps')

snp_near = (snps.Probe_SNPs != '')
snp_near.value_counts()

probe_idx = df_hiv.index.diff(ti(ff))

"""
Explanation: Set up Probe Filters
End of explanation
"""


