<jupyter_start><jupyter_text>This notebook is the analysis supporting the paper titled  _Assessing the impact of general practice software on prescribing safety recommendations for diltiazem, tacrolimus and ciclosporin. A cohort study in English NHS primary care_<jupyter_code>##import libraries neededd
import pandas as pd
import numpy as np
from ebmdatalab import bq, maps, charts
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec

## ensuring the format is consistent for pounds and pence
pd.set_option('display.float_format', lambda x: '%.2f' % x)<jupyter_output><empty_output><jupyter_text>Part 1: Breaches in Guidance by practice/CCG  1a Ciclosporin and Tacrolimus<jupyter_code>sql='''
SELECT
  p.month,
  practice_id,
  p.pct_id,
  numerator,
  denominator,
  calc_value,
  percentile
FROM
  `ebmdatalab.measures.practice_data_ciclosporin` AS p
INNER JOIN
  ebmdatalab.hscic.practices AS prac
ON
  p.practice_id = prac.code
JOIN
  hscic.practice_statistics
ON
  practice_statistics.practice = prac.code
WHERE
  (setting = 4 
  AND denominator > 0)
GROUP BY
p.month,
practice_id,
p.pct_id,
numerator,
denominator,
calc_value,
percentile
'''

df1a = bq.cached_read(sql, csv_path='df1a.csv', use_cache=True)
df1a['month'] = pd.to_datetime(df1a['month'])
df1a.rename(columns={'pct_id':'pct'}, inplace=True) ##prep for maps library whe  pct must be colum name
df1a.head(5)

df1a2018 = df1a.loc[(df1a["month"]>= "2018-01-01") & (df1a["month"]<="2018-12-01")]
df1a2018.head(5)
df1a2018.describe(percentiles=(.1,.90))
df1a2018.sum(numeric_only=True)
#number of practices in the dataset who prescribed at least one Rx 
df1a2018.practice_id.nunique()
#now lets restrict to values where numerator is greater than 0 i.e. some generic Rxs
breachciclo = df1a2018[df1a2018['numerator']>0]
breachciclo.head()
#the number of practices that breached the ciclosporin/tacrolimus in 2018 
breachciclo.practice_id.nunique()
##calculate proprotion that breached
breachciclo.practice_id.nunique()/df1a2018.practice_id.nunique()*100
###here we group for the map
by_pct_df1a = df1a[df1a['month'] == "2019-05-01"].groupby('pct').sum().reset_index()
by_pct_df1a['calc_value'] = by_pct_df1a['numerator'] / by_pct_df1a['denominator'] ##reclaculate the "calc value" as the groupby summed it
by_pct_df1a.head(5)<jupyter_output><empty_output><jupyter_text>1b Diltiazem<jupyter_code>sql='''
SELECT
  p.month,
  practice_id,
  p.pct_id,
  numerator,
  denominator,
  calc_value,
  percentile
FROM
  ebmdatalab.measures.practice_data_diltiazem AS p
INNER JOIN
  ebmdatalab.hscic.practices AS prac
ON
  p.practice_id = prac.code
JOIN
  hscic.practice_statistics
ON
  practice_statistics.practice = prac.code
WHERE
  (setting = 4  ##need to add prac = 4 here
    AND denominator > 0)
GROUP BY
  p.month,
  practice_id,
  p.pct_id,
  numerator,
  denominator,
  calc_value,
  percentile
'''

df1b = bq.cached_read(sql, csv_path='diltiazem.csv', use_cache=True)
df1b['month'] = pd.to_datetime(df1b['month'])
df1b.rename(columns={'pct_id':'pct'}, inplace=True) ##prep for maps library whe  pct must be colum name
df1b.head(5)
df1b2018 = df1b.loc[(df1b["month"]>= "2018-01-01") & (df1b["month"]<="2018-12-01")]
df1b2018.head(5)
df1b2018.describe(percentiles=(.1,.90))
df1b2018.sum(numeric_only=True)
#number of practices in the dataset who prescribed at least one Rx and have a list size >1000?
df1b2018.practice_id.nunique()
#now lets restrict to values where numerator is greater than 0 i.e. some generic Rxs
breachdilt = df1b2018[df1b2018['numerator']>0]
breachdilt.head()
#the number of practices that breached the diltiazem in 2018 
breachdilt.practice_id.nunique()
#calulate proportion of practices that breached
breachdilt.practice_id.nunique()/df1b2018.practice_id.nunique()*100
### here we calculate totals on both measures
total = df1a2018.sum(numeric_only=True) +df1b2018.sum(numeric_only=True)
total
## calculate proprotion of scrips in numerator from total prescritpions
(total['numerator'] / total['denominator'])*100
###here we group for the map
by_pct_df1b = df1b[df1b['month'] == "2019-05-01"].groupby('pct').sum().reset_index()
by_pct_df1b['calc_value'] = by_pct_df1b['numerator'] / by_pct_df1b['denominator'] ##reclaculate the "calc value" as the groupby summed it
by_pct_df1b.head(5)<jupyter_output><empty_output><jupyter_text>Part 2: EHR Work<jupyter_code>##Importing GP EHR Deployment Data
###this will need cleaning up? Seb/ALex?
vendors = pd.read_csv('vendors.csv')
# Clean up the input data
vendors['Principal Supplier'] = vendors['Principal Supplier'].str.strip()
vendors.loc[vendors['Principal Supplier'] == 'INPS', 'Principal Supplier'] = 'Vision'  # seems they changed in 2017
vendors = vendors.loc[vendors['Date'] > '2016-02-01']  # there is some dirty data ("Unknowns") before this

start = pd.to_datetime('2016-01-01')
mid = pd.to_datetime('2017-01-01')
end = pd.to_datetime('2018-12-01')<jupyter_output><empty_output><jupyter_text>2a ciclosporin and tacrolimus measure data linked to EHR provider<jupyter_code>import importlib
importlib.reload(bq)
measures = ['ciclosporin']

def get_data(measure_id):
    sql = """
SELECT
  TRIM(Principal_Supplier) AS supplier,
  m.practice_id,
  m.pct_id,
  Date as month,
  numerator,
  denominator,
  '{measure_id}' AS measure_id
FROM
  measures.practice_data_{measure_id} as m
JOIN
  ebmdatalab.alex.vendors as software #this is where the up to date vendors table is held
ON
  software.ODS = practice_id
  AND DATE(Date) = m.month
JOIN
  hscic.practices
ON
  practices.code = software.ODS
JOIN
  hscic.practice_statistics
ON
  practice_statistics.practice = practices.code
  AND Date = practice_statistics.month
WHERE
total_list_size > 100
  AND practices.status_code = 'A'
  AND denominator > 0
  AND m.month <= "2019-05-01"
ORDER BY
  m.month""".format(measure_id=measure_id)
    import pandas as pd
    df = bq.cached_read(sql, csv_path="ciclosporintacro.csv".format(measure_id))
    return df


df2a = get_data('ciclosporin')
df2a = df2a.replace(['TPP'],'SystmOne')
df2a = df2a.replace(['INPS'],'Vision')
df2a = df2a[df2a.supplier != "Unknown"]
df2a.rename(columns={'pct_id':'pct'}, inplace=True) ##prep for maps library whe  pct must be colum name
df2a.head(5)
df2a['calc_value'] = df2a['numerator'] / df2a['denominator']
df2a['month'] = pd.to_datetime(df2a['month'])
ax = df2a.groupby(['month', 'supplier']).mean()['calc_value'].unstack().plot.line()
plt.ylim(0, 0.25)
plt.legend(loc=(1.04,0))
plt.title(" (a) Total items of generic ciclosporin and tacrolimus preparations  \n as a proportion of total items of all ciclosporin and tacrolimus items by EHR"),
ax.set_ylabel("proportion")
###here we reformat graph to express in term of percentage
ax = df2a.groupby(['month', 'supplier']).mean()['calc_value']
ax = ax*100
ax = ax.unstack().plot.line()
plt.ylim(0, 28)
plt.legend(loc=(1.04,0))
plt.title(" (a) Total items of generic ciclosporin and tacrolimus preparations  \n as a percentage of total items of all ciclosporin and tacrolimus items by EHR"),
ax.set_ylabel("%")
dfccg = df2a.groupby(['month','pct']).sum().reset_index()
dfccg.head(5)<jupyter_output><empty_output><jupyter_text>2b - Diltiazem linked with EHR<jupyter_code>importlib.reload(bq)
measures = ['diltiazem']

def get_data(measure_id):
    sql = """
SELECT
  TRIM(Principal_Supplier) AS supplier,
  m.practice_id,
  m.pct_id,
  Date as month,
  numerator,
  denominator,
  '{measure_id}' AS measure_id
FROM
  measures.practice_data_{measure_id} as m
JOIN
  ebmdatalab.alex.vendors as software #this is where the up to date vendors table is held
ON
  software.ODS = practice_id
  AND DATE(Date) = m.month
JOIN
  hscic.practices
ON
  practices.code = software.ODS
JOIN
  hscic.practice_statistics
ON
  practice_statistics.practice = practices.code
  AND Date = practice_statistics.month
WHERE
total_list_size > 100
  AND practices.status_code = 'A'
  AND denominator > 0
  AND m.month <= "2019-05-01"
ORDER BY
  m.month""".format(measure_id=measure_id)
    import pandas as pd
    df2 = bq.cached_read(sql, csv_path="diltiazem.csv".format(measure_id))
    return df2


df2b = get_data('diltiazem')
df2b = df2b.replace(['TPP'],'SystmOne')
df2b = df2b.replace(['INPS'],'Vision')
df2b = df2b[df2b.supplier != "Unknown"]
df2b.rename(columns={'pct_id':'pct'}, inplace=True) ##prep for maps library whe  pct must be colum name
df2b.head(5)
df2b['calc_value'] = df2b['numerator'] / df2b['denominator']
df2b['month'] = pd.to_datetime(df2b['month'])
df2b.month.nunique
ax1 = df2b.groupby(['month', 'supplier']).mean()['calc_value'].unstack().plot.line()
plt.legend(loc=(1.04,0))
plt.title(" (b) items of generic diltiazem modified-release preparations \n as a proportion of total items of all diltiazem modified-release items (>60mg)")
ax1.set_ylabel("proportion")
## here we reformat the graph to express in terms of percentage
ax1 = df2b.groupby(['month', 'supplier']).mean()['calc_value']
ax1 = ax1*100
ax1 = ax1.unstack().plot.line()
plt.ylim(0, 28)
plt.legend(loc=(1.04,0))
plt.title(" (b) items of generic diltiazem modified-release preparations \n as a percentage of total items of all diltiazem modified-release items (>60mg)")
ax1.set_ylabel("%")
df2b.rename(columns={'pct_id':'pct'}, inplace=True)  # The CCG column must be named 'pct' for the maps function
by_supplier_and_pct = df2b.groupby(['supplier', 'pct']).sum().reset_index()
by_supplier_and_pct['calc_value'] = by_supplier_and_pct['numerator'] / by_supplier_and_pct['denominator']
epr_dec18 = df2b.loc[df2b['month']=='2018-12-01',['practice_id','supplier']]
epr_dec18.head(5)
## here we asess EHR deployment in 2018
ehr_deployment = epr_dec18.groupby('supplier').count()
ehr_deployment.head()<jupyter_output><empty_output><jupyter_text>Ghost Branded Generics can cause prescribing of brands to be "hidden" in prescribing data. We have written about it [extensively on our blog](https://ebmdatalab.net/ghost-branded-generics-why-does-the-cost-of-generic-atorvastatin-vary%EF%BB%BF/). Here we will check how many prescriptions could have been affected.<jupyter_code>sql= '''
SELECT
  DATE(month) AS month,
  pct,
  practice,
  SUM(CASE
      WHEN (possible_savings >= 2 OR possible_savings <=-2) THEN items
    ELSE
    0
  END
    ) AS items
FROM
  ebmdatalab.alex.vw__ghost_generic_measure_paper_2019_07
INNER JOIN
  ebmdatalab.hscic.practices
ON
  practice = code
  AND setting = 4
WHERE
bnf_code LIKE '0802020G0%' OR ##Ciclosporin (Systemic) 
bnf_code LIKE '0802020T0%'  ##Tacrolimus (Systemic) 
GROUP BY
  month,
  pct,
  practice'''

df1a_gbg = bq.cached_read(sql, csv_path='gbg.tacro.ciclo.zip',use_cache=True)
df1a_gbg.head(5)
df1a_gbg.sum(numeric_only = True)
df1a_gbg.loc[(df1a_gbg["month"]>= "2018-01-01") & (df1a_gbg["month"]<="2018-12-01")].sum(numeric_only = True)
sql= '''
SELECT
  DATE(month) AS month,
  pct,
  practice,
  SUM(CASE
      WHEN (possible_savings >= 2 OR possible_savings <=-2) THEN items
    ELSE
    0
  END
    ) AS items
FROM
  ebmdatalab.alex.vw__ghost_generic_measure_paper_2019_07
INNER JOIN
  ebmdatalab.hscic.practices
ON
  practice = code
  AND setting = 4
WHERE
  bnf_code LIKE '0206020C0%AB' OR --Diltiazem HCl_Tab 120mg C/R
  bnf_code LIKE '0206020C0%AC' OR --Diltiazem HCl_Tab 90mg M/R 
  bnf_code LIKE '0206020C0%AD' OR --Diltiazem HCl_Tab 120mg S/R 
  bnf_code LIKE '0206020C0%AE' OR --Diltiazem HCl_Cap 300mg M/R 
  bnf_code LIKE '0206020C0%AF' OR --Diltiazem HCl_Cap 90mg C/R 
  bnf_code LIKE '0206020C0%AG' OR --Diltiazem HCl_Cap 180mg C/R 
  bnf_code LIKE '0206020C0%AH' OR --Diltiazem HCl_Cap 300mg C/R 
  bnf_code LIKE '0206020C0%AI' OR --Diltiazem HCl_Cap 120mg C/R 
  bnf_code LIKE '0206020C0%AK' OR --Diltiazem HCl_Cap 90mg S/R 
  bnf_code LIKE '0206020C0%AL' OR --Diltiazem HCl_Cap 120mg S/R 
  bnf_code LIKE '0206020C0%AM' OR --Diltiazem HCl_Cap 240mg C/R 
  bnf_code LIKE '0206020C0%AN' OR --Diltiazem HCl_Cap 180mg S/R 
  bnf_code LIKE '0206020C0%AP' OR --Diltiazem HCl_Cap 240mg S/R 
  bnf_code LIKE '0206020C0%AS' OR --Diltiazem HCl_Tab 120mg M/R 
  bnf_code LIKE '0206020C0%AT' OR --Diltiazem HCl_Cap 90mg M/R 
  bnf_code LIKE '0206020C0%AU' OR --Diltiazem HCl_Cap 120mg M/R 
  bnf_code LIKE '0206020C0%AV' OR --Diltiazem HCl_Cap 180mg M/R 
  bnf_code LIKE '0206020C0%AW' OR --Diltiazem HCl_Cap 240mg M/R 
  bnf_code LIKE '0206020C0%AX' OR --Diltiazem HCl_Cap 200mg M/R 
  bnf_code LIKE '0206020C0%AY' OR --Diltiazem HCl_Cap 240mg M/R 
  bnf_code LIKE '0206020C0%BA' --Diltiazem HCl_Cap 360mg M/R 
GROUP BY
  month,
  pct,
  practice'''

df1b_gbg = bq.cached_read(sql, csv_path='gbg.diltiazem.zip',use_cache=True)
df1b_gbg.head(5)
df1b_gbg.sum(numeric_only=True)
df1b_gbg.loc[(df1b_gbg["month"]>= "2018-01-01") & (df1b_gbg["month"]<="2018-12-01")].sum(numeric_only = True)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12615649992159322, 'fraction_numerical': 0.05449270816998589, 'mean_word_length': 5.3426156141223275, 'pattern_counts': {'":': 0, '<': 37, '<?xml version=': 0, '>': 47, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '581656', 'license': 'MIT', 'n_tokens_mistral': 5563, 'n_tokens_neox': 4973, 'n_words': 1159, 'path': 'Diltiazem.Tacro.Ciclosporin.Paper.ipynb', 'repo_name': 'HDRUK/software_differences_paper', 'size': '12754'}