Diabetes Reports | Python¶

# define diabetes patients Ids
diab_patient_ids = conditions[conditions['CODE'] == 44054006]['PATIENT'].unique()
print("Number of diabetes patients:", diab_patient_ids.size)

diab_patients = patients[patients['Id'].isin(diab_patient_ids)]
print(diab_patients.shape)
diab_patients.head(2)

# extract conditions for diabetes patients
diab_conditions = conditions[conditions['PATIENT'].isin(diab_patient_ids)]
diab_conditions.head(2)

# for inpatients
inpatient_ids = encounters[encounters['ENCOUNTERCLASS'] == 'inpatient']['PATIENT'].unique()
diab_inpatient_ids = np.intersect1d(diab_patient_ids, inpatient_ids)
print("number of inpatient patients:", diab_inpatient_ids.size)
diab_inpatient_patients = patients[patients['Id'].isin(diab_inpatient_ids)]

# for diabetes patient deaths
death_ids = patients[patients.DEATHDATE.notna()].Id.unique()
diab_death_ids = np.intersect1d(diab_patient_ids, death_ids)
print("number of diabetes patient deaths:", diab_death_ids.size)

import re
# Make an explicit copy if needed
diab_conditions = diab_conditions.copy()

# Normalize to lowercase using .loc
diab_conditions.loc[:, 'DESCRIPTION'] = diab_conditions['DESCRIPTION'].str.lower()

# Define diabetes-related keywords
diabetes_related_keywords = [
    'diabetes',
    'prediabetes',
    'hyperglycemia',
    'hypoglycemia',
    'neuropathy',
    'nephropathy',
    'retinopathy',
    'vision',
    'foot ulcer',
    'amputation',
    'cardiovascular',
    'metabolic syndrome',
    'hypertriglyceridemia',
    'obesity',
    'hypertension',
    'anemia',
    'infection',
    'slow healing',
    'renal failure',
    'fatigue',
    'blurred vision'
]

# Build regex pattern with word boundaries
pattern = '|'.join([r'\b' + re.escape(word) + r'\b' for word in diabetes_related_keywords])

# Create mask with regex
mask = diab_conditions['DESCRIPTION'].str.contains(pattern, regex=True)
diabetes_outcomes = diab_conditions[mask]

# Further filter rows where 'DESCRIPTION' contains '(disorder)'
disorder_mask = diabetes_outcomes['DESCRIPTION'].str.contains(r'\(disorder\)', na=False)
diabetes_outcomes = diabetes_outcomes[disorder_mask]
diabetes_outcomes['DESCRIPTION'] = diabetes_outcomes['DESCRIPTION'].str.replace(
    r'\(disorder\)', '', regex=True).str.strip()

# Calculate value counts and percentages
counts = diabetes_outcomes['DESCRIPTION'].value_counts().head(12)
percentages = round((counts / diab_patient_ids.size) * 100, 2)

# Combine counts and percentages into a DataFrame
diabetes_outcomes = pd.DataFrame({
    'DESCRIPTION' : counts.index,
    'Count': counts.values,
    'Percentage': percentages.values
    })

# Display the DataFrame
print(diabetes_outcomes)

counts = diabetes_outcomes['Count']

plt.figure(figsize=(9, 5))
ax = counts.plot(kind='barh', color='teal')
plt.xlabel("Frequency")
plt.title("Top 10 Diabetes-Related Outcomes")
plt.gca().invert_yaxis()
plt.yticks(ticks=range(len(diabetes_outcomes)), labels=diabetes_outcomes['DESCRIPTION'])

# Add value labels inside the bars
for index, value in enumerate(counts):
    plt.text(value, index, str(value), va='center', ha='left', fontsize=10)

plt.tight_layout()
plt.show()

diab_conditions = conditions[conditions['CODE'] == 44054006].copy()
diab_conditions.rename(columns = {'START' : 'diab_date'}, inplace=True)
diab_conditions.head(2)

retino_conditions = diab_conditions.merge(conditions, on='PATIENT')
retino_conditions.head(2)

# Filter for diabetic retinopathy descriptions
retino_conditions = retino_conditions[(retino_conditions['CODE_y'] == 422034002) | 
                                      (retino_conditions['CODE_y'] == 1551000119108)]
# Ensure START and diab_date columns are in datetime format
retino_conditions['START'] = pd.to_datetime(retino_conditions['START'], errors='coerce')
retino_conditions['diab_date'] = pd.to_datetime(retino_conditions['diab_date'], errors='coerce')

# Calculate diab_years
retino_conditions['diab_years'] = (((retino_conditions['START'] - retino_conditions['diab_date']).dt.days) / 365.25).round(2)
retino_conditions.head(2)

# Define age group bins and labels
bins = [0, 10, 20, 30, 40]
labels = ['0-9', '10-19', '20-29', '30+']

# Create age groups
retino_conditions['year_group'] = pd.cut(retino_conditions['diab_years'], bins=bins, labels=labels, right=False)
retino_conditions.pivot_table(
    values='DESCRIPTION_y',
    index='year_group',
    columns='CODE_y',
    aggfunc='count',
    observed=False
)

observations['DESCRIPTION'][observations['PATIENT'].isin(diab_inpatient_ids)].value_counts().head()

labs_code = ["2339-0", "4548-4", "2571-8", "38483-4", "2947-0", "6298-4"]
# Filter observations for the specified lab codes and inpatient IDs
diab_observations = observations[
    observations['CODE'].isin(labs_code) & 
    observations['PATIENT'].isin(diab_inpatient_ids)
]
diab_observations.head(2)

diab_patients_labs = diab_conditions.merge(diab_observations, on='PATIENT')
diab_patients_labs['DATE'] = pd.to_datetime(diab_patients_labs.DATE)
diab_patients_labs['labs_days'] = (diab_patients_labs.DATE - diab_patients_labs.diab_date).dt.days
diab_patients_labs['years'] = round((diab_patients_labs['labs_days'] / 365.25), 2)
diab_patients_labs['VALUE'] = pd.to_numeric(diab_patients_labs['VALUE'], errors='coerce')
diab_patients_labs = diab_patients_labs[diab_patients_labs['years'] > 0]
diab_patients_labs.rename(columns={'CODE_y': 'CODE'}, inplace=True)
diab_patients_labs.rename(columns={'DESCRIPTION_y': 'DESCRIPTION'}, inplace=True)
diab_patients_labs = diab_patients_labs[['diab_date', 'PATIENT', 'CODE', 'DESCRIPTION', 'VALUE', 'years']]
diab_patients_labs['hospitalized'] = diab_patients_labs.PATIENT.isin(inpatient_ids)
diab_patients_labs.head(2)

loinc_to_display = {'CODE = 2339-0': 'Glucose', 'CODE = 4548-4': 'Hemoglobin A1c/Hemoglobin.total in Blood', 'CODE = 2571-8': 'Triglycerides',  'CODE = 38483-4': 'Creatinine', 'CODE = 2947-0': 'Sodium', 'CODE = 6298-4': 'Potassium'}

catplt = sns.catplot(x="years", y="VALUE", kind="box", col='CODE', 
            col_wrap=2, sharey=False, sharex=False, data = diab_patients_labs[diab_patients_labs['years']<= 3])

# Set the title for the entire figure
catplt.fig.suptitle("Labs Values of Diabetes Patients", y=1.02)

for axis in catplt.fig.axes:
    axis.xaxis.set_major_formatter(ticker.FormatStrFormatter('%d'))
    axis.xaxis.set_major_locator(ticker.MultipleLocator(base=4))
    axis.set_title(loinc_to_display[axis.get_title()])
        
plt.show()

# --- Step 1: Prepare the A1C Data ---
diab_A1C = diab_observations[diab_observations['CODE'] == '4548-4'].copy()
diab_A1C['DATE'] = pd.to_datetime(diab_A1C['DATE'], errors='coerce')
diab_A1C['A1C'] = pd.to_numeric(diab_A1C['VALUE'], errors='coerce')
diab_A1C = diab_A1C[['DATE', 'PATIENT', 'A1C']].dropna()

# --- Step 2: Find Baseline Date per Patient ---
baseline_dates = diab_A1C.groupby('PATIENT')['DATE'].min().reset_index()
baseline_dates.rename(columns={'DATE': 'BASELINE_DATE'}, inplace=True)

# --- Step 3: Merge and Calculate Days Since Baseline ---
diab_A1C = diab_A1C.merge(baseline_dates, on='PATIENT')
diab_A1C['DAYS_SINCE_BASELINE'] = (diab_A1C['DATE'] - diab_A1C['BASELINE_DATE']).dt.days

# --- Step 4: Assign Timepoints ---
def assign_timepoint(days):
    if days == 0:
        return 'baseline'
    elif 80 <= days <= 100:
        return '3_month'
    elif 160 <= days <= 200:
        return '6_month'
    else:
        return None

diab_A1C['TIMEPOINT'] = diab_A1C['DAYS_SINCE_BASELINE'].apply(assign_timepoint)
diab_A1C = diab_A1C.dropna(subset=['TIMEPOINT'])

# --- Step 5: Aggregate Statistics ---
summary = diab_A1C.groupby('TIMEPOINT')['A1C'].agg(
    N='count',
    mean='mean',
    sd='std'
).reset_index()

# --- Output Results ---
order = ['baseline', '3_month', '6_month']
summary = summary.set_index('TIMEPOINT').loc[order].reset_index()
print(summary)

plt.figure(figsize=(10, 6))
sns.boxplot(x='TIMEPOINT', y='A1C', data=diab_A1C, order=['baseline', '3_month', '6_month'])
plt.xticks(rotation=0, ha='right')
plt.title('A1C at Time Points for Diabetes Patients')
plt.xlabel('Time Points')
plt.ylabel('A1C')
plt.tight_layout()
plt.show()

from scipy.stats import f_oneway

# Filter only valid timepoints with enough data
valid_timepoints = diab_A1C['TIMEPOINT'].unique()
groups = [diab_A1C.loc[diab_A1C['TIMEPOINT'] == tp, 'A1C'] for tp in valid_timepoints if diab_A1C.loc[diab_A1C['TIMEPOINT'] == tp, 'A1C'].size > 1]

# Perform one-way ANOVA
anova_result = f_oneway(*groups)
print("ANOVA F-statistic:", anova_result.statistic)
print("ANOVA p-value:", anova_result.pvalue)

diab_care = care_plans[care_plans['PATIENT'].isin(diab_patient_ids)]
# Normalize REASONDESCRIPTION to lowercase
diab_care.loc[:, 'REASONDESCRIPTION'] = diab_care['REASONDESCRIPTION'].str.lower()
# Use the existing diabetes-related keywords pattern
reason_mask = diab_care['REASONDESCRIPTION'].str.contains(pattern, regex=True, na=False)

# Filter the medications DataFrame
diabetes_related_reasons = diab_care[reason_mask]

# Calculate value counts and percentages
counts = diabetes_related_reasons['REASONDESCRIPTION'].value_counts()
percentages = round((counts / diab_patient_ids.size) * 100, 2)

# Combine counts and percentages into a DataFrame
diabetes_careplans = pd.DataFrame({
    'DESCRIPTION' : counts.index,
    'Count': counts.values,
    'Percentage': percentages.values
    })
diabetes_careplans

diab_conditions = conditions[conditions['CODE'] == 44054006].copy()
diab_conditions.rename(columns = {'START' : 'diab_date'}, inplace=True)
diab_conditions.head(2)

dp = diab_conditions.merge(patients, how='left', left_on='PATIENT', right_on='Id')
dp = dp[['diab_date', 'PATIENT', 'BIRTHDATE', 'DEATHDATE', 'GENDER']]
dp['diab_date'] = pd.to_datetime(dp['diab_date'])
dp['BIRTHDATE'] = pd.to_datetime(dp['BIRTHDATE'])
dp['DEATHDATE'] = pd.to_datetime(dp['DEATHDATE'], errors='coerce')  # Handle NaT for missing dates
dp['admit'] = dp.PATIENT.isin(inpatient_ids)
dp['death'] = dp.DEATHDATE.notna()
dp['death_age'] = ((dp['DEATHDATE'] - dp['BIRTHDATE']).dt.days / 365.25).round(2)  
dp['diab_age'] = ((dp['diab_date'] - dp['BIRTHDATE']).dt.days / 365.25).round(2)
print("mean age at death:", dp['death_age'].mean().round(2), "years-old")
dp.head(3)

pd.crosstab(dp['GENDER'], dp['death'])

# Create a cross-tabulation with percentages
cross_tab = pd.crosstab(dp['GENDER'], dp['death'], normalize='index') * 100

# Display the cross-tabulation
print(cross_tab.round(2))

dp_death = dp[dp['DEATHDATE'].notna()]
dp_death.head(2)

death = pd.DataFrame(dp_death['death'].value_counts()).reset_index()
death.columns = ['death', 'count']
deaths = death.loc[death['death'] == True, 'count'].iloc[0] # counting the deaths
death_rate = (deaths/diab_patient_ids.size)
print("Diabetes death rate is:", round(death_rate*100, 2), "%")

# Filter only the rows where death occurred
diabetes_deaths = dp[dp['death'] == 1]

# Plot the histogram
plt.figure(figsize=(8, 4))
plt.hist(diabetes_deaths['diab_age'], bins=20, color='skyblue', edgecolor='black')
plt.xlabel('Age at Diabetes Diagnosis (Years)')
plt.ylabel('Number of Deaths')
plt.title('Number of Diabetes Deaths by Age at Diagnosis')
plt.grid(axis='y', linestyle='--', alpha=0.7)
plt.tight_layout()
plt.show()

dp.pivot_table('death', index='admit', columns = 'GENDER').round(3)

diab_meds = medications[medications['PATIENT'].isin(diab_patient_ids)]
count = diab_meds['DESCRIPTION'].value_counts()
percentages = round((count / count.sum()) * 100, 2)

# Combine counts and percentages into a DataFrame
diabetes_meds = pd.DataFrame({
    'Medications' : count.index,
    'Count': count.values,
    'Percentage': percentages.values
    })
diabetes_meds.head(14)

# use function size()
diab_meds.groupby('DESCRIPTION').size().reset_index(name='count').sort_values('count', ascending=False).head()

diab_encounters = encounters[encounters['ENCOUNTERCLASS'] == 'inpatient']
diab_encounters = diab_encounters.copy()
diab_encounters['stay_days'] = (
    pd.to_datetime(diab_encounters['STOP']) - pd.to_datetime(diab_encounters['START'])
) / np.timedelta64(1, 'D')

diab_encounters['stay_days'].describe()

# Remove rows with the abnormal maximum stay_days value
diab_encounters = diab_encounters[diab_encounters['stay_days'] != diab_encounters['stay_days'].max()]
(diab_encounters['stay_days'].sum() / diab_inpatient_ids.size).round(2)