Pharmacy In Practice

David Fong

4/Oct/2020

Abstract

Pharmacists in general practice are expected to improve medication efficacy and reduce medication complications and adverse effects. In order for pharmacists in practice to be sustainable, practice revenue needs to be taken into consideration. In this short study, conducted during a period of considerable practice and societal upheaval (December 2019 to August 2020) due to a coronavirus outbreak, a pharmacist in practice was able to reduce the overall medication count (-0.52 per patient) and increase the billings (‘standard’ +0.86, ‘nursing’ +0.29, ‘case conference’ +0.03 and ‘medication review’ +0.09) of patients the pharmacist contacted.

Patients who were contacted by the pharmacists were more likely to have frequent contacts with the practice, have more medications in their medication list and more likely to have been admitted to the hospital prior to the study period. These characteristics are associated with high risk of future admission to hospital.

Additional revenue per unique patient contact is estimated at $90.14. One hundred and seventy unique patients were contacted, resulting in a total estimated revenue increase of $15323.80. Financial viability will require measures to ensure the pharmacist has an adequate number of patient contacts and measures to ensure follow-through in claimable items.

1 Introduction

General practice involves a significant, and ever increasing, number of medical therapeutics, often prescribed to patients with increasing age and co-morbidities. At the coHealth Kensington site, many patients also come from a background with reduced English communication ability and overall health literacy. Many have been displaced from another country and have moved to Australia at a mature age. coHealth Kensington is located close to areas with a high density of public housing and sees patients who are homeless or at risk of homelessness. The site also provides mental health support. A high proportion of patients have mental health conditions, including low-frequency disorders such as schizophrenia.

This is on the background of a local population which is ageing along with the general Australian community.

The pharmacy in practice program, provided by the Primary Health Network, supported the placement of an ‘in-house’ pharmacist at the Kensington site for several days a week. The pharmacist started at the Kensington site on 2nd December 2019.

A pharmacist in practice is expected to improve medication efficacy (through appropriate use, simplification and patient education) and reduce medication complications (through simplification, reduced misunderstanding and review of potential medication interactions and contra-indications). For a pharmacist in practice to be sustainable, practice revenue needs to be taken into consideration.

This evaluation examines whether contacts by the pharmacist with patients reduced overall medication counts (as listed in the electronic medical record) and contributed to the billings/finances of the practice.

The evaluation is complicated by the effects of a coronavirus outbreak throughout much of 2020, which severely curtailed the number of patients who physically came to the clinic site and overall changed the way the practice operated.

library(readxl)    # read Excel spreadsheets (with gritted teeth...)
library(summarytools)
library(magrittr)
library(dplyr)     # data pliers
library(tidyr)
library(lubridate)
library(ggplot2)   # plotting
library(cowplot)   # side-by-side plot
library(lme4)      # glmer
library(microsynth) # synthetic models
library(sjPlot)    # regression table display
library(kableExtra) # table display

2 Sources of measurement

The active patient list
- generated by dMeasure/GPstat!, see Appendix C
- Patients billed on three or more separate days between 2017-10-01 and 2019-09-30 inclusive
Patient details (date of birth, age calculation)
- generated by PenCS CAT4
The dated list of billings for patient in the ‘active’ list
- generated by dMeasure/GPstat!
List of telephone and ‘in-person (in surgery)’ contact by the pharmacist with patients
- generated by dMeasure/GPstat!
Medication counts of Kensington pateints 2019-12-01 and 2020-08-01
- generated by PenCS CAT4

### read in data

active_patient <- read.csv("active_patient.csv")
patient_details <- read_excel("Kensington20200901_all.xlsx", skip = 3)
active_patient_billings <- read.csv("active_patient_billings.csv")
pharmacist_contacts <- read.csv("LauraDeanVisits.csv")
discharge <- read.csv("KensingtonDischarge_active.csv")


# change various numericals into 'categoricals'
active_patient <- active_patient %>%
  mutate(
    InternalID = factor(InternalID)
  ) %>%
  left_join(
    patient_details %>%
      # create DOB, AgeYears (age in years as of 2/Dec/2019)
      # and AgeGroup (20 year age groups)
      select(ID, `D.O.B (Age)`) %>%
      rename(InternalID = ID) %>%
      mutate(
        DOB = as.Date(substr(`D.O.B (Age)`, 1, 10), format = "%d/%m/%Y"),
        AgeYears = as.numeric(floor(
          time_length(interval(DOB, as.Date("2019/12/02")), "years")
        ))
      ) %>%
      mutate(AgeGroup = pmax(floor(AgeYears/20), 0)) %>%
      select(InternalID, DOB, AgeYears, AgeGroup),
    by = "InternalID"
  ) %>%
  filter(!is.na(DOB)) 
# patients without DOB as of 1/Sep/2020 are actually *not* active as of
# that timepoint e.g. passed away earlier
      
# also change some date 'strings' to true 'dates'
active_patient_billings <- active_patient_billings %>%
  mutate(
    InternalID = factor(InternalID),
    MBSItem = factor(MBSItem),
    ServiceDate = as.Date(ServiceDate)
  ) %>%
  filter(Paid > 0) # only those billings which have actually been paid
# these are not truly numeric values, they are actually 'categoricals'

pharmacist_contacts <- pharmacist_contacts %>%
  mutate(
    InternalID = factor(InternalID),
    VisitType = factor(VisitType),
    VisitDate = as.Date(VisitDate)
  )

discharge <- discharge %>%
  mutate(
    InternalID = factor(InternalID),
    CorrespondenceDate = as.Date(CorrespondenceDate)
  )

medication_count <- list()
medication_count[[1]] <- read_excel(
  "PenCSMedicationCount20191201.xlsx", sheet = "ReidentifyReport",
  skip = 3) %>%
  mutate(DOB = as.Date(substr(`D.O.B (Age)`, 1, 10), format = "%d/%m/%Y"),
         Sex = factor(Sex)) %>%
  # age was not a generally useful predictor on testing
  #  it does help predict CDM and HA, but had no effect on results
  select(ID, `Med. Count`, Sex) %>%
  rename(InternalID = ID, MedCount = `Med. Count`) %>%
  mutate(InternalID = factor(InternalID), MedCount = as.numeric(MedCount))
medication_count[[2]] <- read_excel(
  "PenCSMedicationCount20200801.xlsx", sheet = "ReidentifyReport",
  skip = 3) %>%
  select(ID, `Med. Count`) %>%
  rename(InternalID = ID, MedCount = `Med. Count`) %>%
  mutate(InternalID = factor(InternalID), MedCount = as.numeric(MedCount))

3 Specification of the measured information (data)

3.1 Patient contacts (and potential contact) with pharmacist

clinicVisits - number of times patient physically came to the clinic
- additional sub-categories into clinicVisit1 and clinicVisit2 - number of times the patient came to the clinic in period 1 and period 2
- Period 1 : the twelve months before 2019-12-02 (2nd December 2019)
- Period 2 : the nine months after 2019-12-2 (2nd December 2019), including 2nd December 2019
pharmacistExposure - number of times patient physically came to the clinic on the day the pharmacist was present.
DirectContact - TRUE/FALSE - pharmacist contacted patient either in surgery or via telephone.

physical_visit_MBSItem <- c("23", "36", "3", "44", "10997", "721", "723", "732",
                            "2713", "2715", "73806", "16500", "707", "705", "2517", "2521",
                            "2525", "703", "30071", "2501", "2504")

pharmacist_dates <- read.csv("LauraDeanDays.csv") %>% pull(Date) %>% as.Date()

active_patient_inClinic <- active_patient_billings %>%
  mutate(physicalVisit = MBSItem %in% physical_visit_MBSItem) %>%
  group_by(InternalID, ServiceDate) %>% # group items on the same day/patient together
  summarise(ClinicVisit = any(physicalVisit)) %>%
  ungroup() %>%
  mutate(pharmacistDay = ClinicVisit & ServiceDate %in% pharmacist_dates)
# pharmacistDay is TRUE is patient came to clinic physically AND
#  pharmacist was present on the same day

# first day of pharmacist visit was 2019-12-02 (2nd December 2019)

ClinicExposure_df <- active_patient_inClinic %>%
  group_by(InternalID) %>%
  summarise(clinicVisits = sum(ClinicVisit),
            clinicVisit1 = sum(ClinicVisit & ServiceDate >= (as.Date("2019-12-02") - months(12)) & ServiceDate < as.Date("2019-12-02")),
            # clinicVisit1 'pre' (12 month) period before 2nd of December, the number of 'physical' visits
            clinicVisit2 = sum(ClinicVisit & ServiceDate >= as.Date("2019-12-02") & ServiceDate <= (as.Date("2019-12-02") + months(9))),
            # clinicVisit2 'post' (9 month) period after and including 2nd December, the number of 'physical' visits
            pharmacistExposure = sum(pharmacistDay)) %>%
  ungroup()

DirectContact_df <- pharmacist_contacts %>%
  group_by(InternalID) %>%
  summarise(DirectContact = any(VisitType %in% c("Surgery", "Telephone"))) %>%
  ungroup()

# add the 'exposures' to the table
#  ClinicExposure - patient came on same day that pharmacist was present
#  DirectContact - pharmacist made contact with patient, either by telephone or in-person

active_patient_exposure <- active_patient %>%
  left_join(ClinicExposure_df, by = "InternalID") %>%
  left_join(DirectContact_df, by = "InternalID") %>%
  replace_na(list(DirectContact = FALSE))

3.2 Medication

Each patient has a MedCount.
MedCount is measured at two time points 1 and 2
- MedCount1 : 1st December 2019
- MedCount2 : 1st August 2020
a MedCountDelta is calculated, subtracting the count in 1 from 2.

active_patient_medication <- active_patient_exposure %>%
  left_join(medication_count[[1]] %>%
              rename(MedCount1 = MedCount),
            by = "InternalID") %>%
  left_join(medication_count[[2]] %>%
              rename(MedCount2 = MedCount),
            by = "InternalID") %>%
  mutate(MedCountDelta = MedCount2 - MedCount1)

3.3 Discharge

discharges - total number of times patient has a discharge notification
- additional sub-categories into discharge1 and discharge2 - number of times the patient came to the clinic in period 1 and period 2
- Period 1 : the twelve months before 2019-12-02 (2nd December 2019)
- Period 2 : the nine months after 2019-12-2 (2nd December 2019), including 2nd December 2019
dataframe active_patient_discharge
- each InternalID patent has rows of weekly row_dates
- pharmacyContact is TRUE if pharmacist contacted before that date
- admission is TRUE is a discharge correspondence received just after that date.

discharge_selection <- discharge %>%
  # just the discharge notices 12 months before or 9 months after 2/Dec/2019
  filter(CorrespondenceDate >= (as.Date("2019-12-02") - months(12)) &
           CorrespondenceDate < as.Date("2019-12-02") + months(9))
  
Discharge_df <- discharge_selection %>%
  group_by(InternalID) %>%
  summarise(discharges = n(),
            discharge1 = sum(CorrespondenceDate >= (as.Date("2019-12-02") - months(12)) & CorrespondenceDate < as.Date("2019-12-02")),
            discharge2 = sum(CorrespondenceDate >= as.Date("2019-12-02") & CorrespondenceDate < (as.Date("2019-12-02") + months(9)))
            ) %>%
  # number of discharges before 2/Dec/2019 and after 2/Dec/2019
  ungroup() %>%
  replace_na(list(discharge = 0, discharge1 = 0, discharge2 = 0))

record_dates <- seq.Date(as.Date("2019-12-02") - months(12),
                         as.Date("2019-12-02") + months(9),
                         by = "week")
# the 'recording dates' for discharges and pharmacy contacts

active_patient_discharge <- active_patient_medication %>%
  left_join(Discharge_df, by = "InternalID") %>%
  select(InternalID, clinicVisit1, clinicVisit2, Sex, MedCount1,
         discharge1, discharge2, AgeYears, AgeGroup) %>%
  replace_na(list(discharge1 = 0, discharge2 = 0,
                  MedCount1 = 0, clinicVisit1 = 0, clinicVisit2 = 0)) %>%
  mutate(row_date = list(record_dates)) %>%
  unnest(row_date) %>% # this creates a row for each of the dates (92 weeks) %>%
  mutate(admission = 0, pharmacistContact = 0)

for (i in rownames(pharmacist_contacts)) {
  if (pharmacist_contacts[i, "VisitType"] %in% c("Surgery", "Telephone")) {
    intID <- as.numeric(as.character(pharmacist_contacts[i, "InternalID"]))
    visitDate <- pharmacist_contacts[i, "VisitDate"]
    active_patient_discharge <- active_patient_discharge %>%
      mutate(pharmacistContact = if_else(
        as.numeric(as.character(InternalID)) == intID & row_date > visitDate,
        1,
        pharmacistContact))
    # marked as pharmacist contact in the week *after* the actual contact date
  }
}

for (i in rownames(discharge_selection)) {
  intID <- as.numeric(as.character(discharge_selection[i, "InternalID"]))
  admission_date <- max(record_dates[record_dates <= discharge_selection[i, "CorrespondenceDate"]])
  # choose a date *less* than the correspondence date
  active_patient_discharge <- active_patient_discharge %>%
    mutate(admission = if_else(
      as.numeric(as.character(InternalID)) == intID & row_date == admission_date,
      1,
      admission
    ))
}

active_patient_discharge <- active_patient_discharge %>% 
  mutate(timeVar = match(row_date, record_dates)) # create ranked version of dates

3.4 Billings numbers and frequencies

The billings numbers, and frequencies, listed in declining order.

summary(active_patient_billings$MBSItem)

##   10990      23      36       3      44   10997   91809     721     732     723 
##   21341    8096    6515    1824    1590    1270    1103     900     828     795 
##   10960   10981    2713   91810   10962   74990    2715      37   73806   16500 
##     586     560     534     454     291     181     172     161     154     141 
##       0    2712   91811     116     707    2717   91795      47      24    2517 
##     127     127     118     116     110     106      96      78      73      67 
##   41647   73811   90001   93201   91800   92068     699     705     133    2521 
##      62      60      59      53      46      46      45      45      43      41 
##   92069     132   90043     197     703     735   11700   93203   30071   91801 
##      41      40      40      34      28      28      27      22      21      21 
##   14206     747     900   30062   90035   92072   92128    2525     193     739 
##      20      19      19      19      18      18      18      17      16      14 
##   11506   91835     750   80010       4     743   92127   92129      43     110 
##      14      13      12      12      11      11      11      11       9       9 
##     758   30219     701    2501   10956   31366   91802   92126    5049   11702 
##       8       8       7       7       7       7       7       7       6       6 
##    6057   11505   31365    2504    6058   35503   90051   91825    2546    5028 
##       5       5       5       4       4       4       4       4       3       3 
##   18242   41500   73805   92024   92028      35     199     903    2507 (Other) 
##       3       3       3       3       3       2       2       2       2      46

3.5 Billing counts

active_patient_billings_period <- list()
active_patient_billings_period[[1]] <- active_patient_billings %>% 
  filter(ServiceDate >= (as.Date("2019-12-02") - months(12)),
         ServiceDate < (as.Date("2019-12-02")))
active_patient_billings_period[[2]] <- active_patient_billings %>%
  filter(ServiceDate >= (as.Date("2019-12-02")),
         ServiceDate < (as.Date("2019-12-02") + months(9)))

Grouped into Standard, CDM, HA, DMMR, Conference, Nurse

Standard - 3, 23, 36, 44, 2713 etc. and telephone/telehealth equivalents
CDM - 721, 723, 732 and telephone/telehealth equivalents
Nurse - 10997 and telephone/telehealth equivalents
HA - 701, 703, 705, 707
DMMR - 900
Conference - 735, 739, 743, 747, 750, 758 case conference
Total - sum of all the above groups

active_patient_billings_counts <- lapply(
  active_patient_billings_period,
  function(x) {
    # count the number of times each billing category occurred in each 9 month period
    x %>%
      group_by(InternalID) %>%
      summarise(Standard = sum(MBSItem %in% c("3", "23", "36", "44", "2713",
                                              "91790", "91800", "91801", "91802", "92115",
                                              "91795", "91809", "91810", "91811", "92127")),
                CDM = sum(MBSItem %in% c("721", "723", "732",
                                         "92024", "92025", "92027",
                                         "92068", "92069", "92071")),
                HA = sum(MBSItem %in% c("701", "703", "705", "707")),
                DMMR = sum(MBSItem %in% c("900")),
                Conference = sum(MBSItem %in% c("735", "739", "743",
                                                "747", "750", "758")),
                Nurse = sum(MBSItem %in% c("10997", "93201", "93203"))) %>%
      ungroup()
  }
)

Each billing category is counted in the two (1 and 2) time periods

Period 1 : the twelve months before 2019-12-02 (2nd December 2019)
Period 2 : the nine months after 2019-12-2 (2nd December 2019), including 2nd December 2019

Each billing category has a Delta count, calculated by subtracting ‘first (1)’ period from the ‘second (2)’ period count. The Delta count can be, and often was, negative because:

the second time period (9 months) is shorter than the first time period (12 months)
the second period included an extensive lock-down period due to a coronavirus pandemic. Resulting in both reduced access and qualitative changes in service access.

active_patient_delta <- active_patient_medication %>%
  left_join(active_patient_billings_counts[[1]] %>%
              rename(Standard1 = Standard,
                     CDM1 = CDM,
                     HA1 = HA,
                     DMMR1 = DMMR,
                     Conference1 = Conference,
                     Nurse1 = Nurse),
            by = "InternalID") %>%
  left_join(active_patient_billings_counts[[2]] %>%
              rename(Standard2 = Standard,
                     CDM2 = CDM,
                     HA2 = HA,
                     DMMR2 = DMMR,
                     Conference2 = Conference,
                     Nurse2 = Nurse),
            by = "InternalID") %>%
  replace_na(list(Standard1 = 0, Standard2 = 0,
                  CDM1 = 0, CDM2 = 0,
                  HA1 = 0, HA2 = 0,
                  DMMR1 = 0, DMMR2 = 0,
                  Conference1 = 0, Conference2 = 0,
                  Nurse1 = 0, Nurse2 = 0)) %>%
  mutate(Total1 = Standard1 + CDM1 + HA1 + DMMR1 + Conference1 + Nurse1,
         Total2 = Standard2 + CDM2 + HA2 + DMMR2 + Conference2 + Nurse2) %>%
  filter(Total1 > 0,
         Total2 > 0) %>%
  # only patients who were billed *at least once* in BOTH period 1 AND period 2
  mutate(StandardDelta = Standard2 - Standard1,
         CDMDelta = CDM2 - CDM1,
         HADelta = HA2 - HA1,
         DMMRDelta = DMMR2 - DMMR1,
         ConferenceDelta = Conference2 - Conference1,
         NurseDelta = Nurse2 - Nurse1,
         TotalDelta = Total2 - Total1)

4 Patient summary

Number of unique patients who have had at least had 3 billed contacts at Kensington over the past 24 months.

(Restricted to those who have had at least one billing in the twelve months prior to 2nd December 2019, and at least one billing in the nine months after 2nd December 2019.

nrow(active_patient_delta)

## [1] 1196

Number of patients who attended the clinic (physically) on the days the pharmacist was present:

nrow(active_patient_delta %>% filter(pharmacistExposure > 0))

## [1] 826

Number of patients directly contacted by pharmacist, either in person or by telephone:

nrow(active_patient_delta %>% filter(DirectContact > 0))

## [1] 168

4.1 Number of physical visits to clinic

Period 1 - 12 month period prior to 2nd December 2019

Period 2 - 9 month period after and including 2nd December 2019

descr(active_patient_delta %>% 
        select(clinicVisit1, clinicVisit2),
      stats = c("mean", "sd", "min", "q1", "med", "q3", "max"),
      style = "rmarkdown", headings = FALSE
)

	clinicVisit1	clinicVisit2
Mean	7.81	3.82
Std.Dev	5.64	3.43
Min	1.00	0.00
Q1	4.00	1.00
Median	6.00	3.00
Q3	10.00	5.00
Max	48.00	36.00

4.2 Hospital admissions and discharges

Hospital admissions and discharges in period 2 (2nd December 2019 to 2nd September 2020) are positively correlated with the number of admissions and discharges in period 1 (the twelve months before 2nd December 2019), the number of clinic visits in period 1, the number of medications in the patient list in period 1 and male sex.

Number of admissions is also related to the ‘age group’. In the model below, age groups are twenty year age groups, starting at age zero years.

discharge_model <- lm(
  discharge2 ~ discharge1 + clinicVisit1 + MedCount1 + factor(AgeGroup) + Sex,
  data = active_patient_discharge
)
tab_model(discharge_model)

	discharge2
Predictors	Estimates	CI	p
(Intercept)	-0.04	-0.05 – -0.03	<0.001
discharge1	0.24	0.24 – 0.24	<0.001
clinicVisit1	0.02	0.01 – 0.02	<0.001
MedCount1	0.01	0.01 – 0.01	<0.001
AgeGroup [1]	0.02	0.01 – 0.03	0.001
AgeGroup [2]	-0.04	-0.05 – -0.03	<0.001
AgeGroup [3]	0.01	-0.00 – 0.02	0.242
AgeGroup [4]	-0.04	-0.06 – -0.03	<0.001
Sex [M]	0.02	0.02 – 0.03	<0.001
Observations	137632
R² / R² adjusted	0.192 / 0.192

4.3 Billings (of select categories)

4.3.1 ‘Total’ billings (of all selected categories)

descr(active_patient_delta %>%
          select(Total1, Total2),
      stats = c("mean", "sd", "min", "q1", "med", "q3", "max"),
      style = "rmarkdown", headings = FALSE
)

	Total1	Total2
Mean	8.68	6.12
Std.Dev	6.44	5.51
Min	1.00	1.00
Q1	4.00	2.00
Median	7.00	4.00
Q3	12.00	8.00
Max	56.00	45.00

4.3.2 ‘Standard’ billings

Item A/B/C/D, including telehealth/telephone.
Period 1 - 12 month period prior to 2nd December 2019
Period 2 - 9 month period after and including 2nd December 2019

descr(active_patient_delta %>%
          select(Standard1, Standard2),
      stats = c("mean", "sd", "min", "q1", "med", "q3", "max"),
      style = "rmarkdown", headings = FALSE
)

	Standard1	Standard2
Mean	7.05	5.00
Std.Dev	5.35	4.42
Min	0.00	0.00
Q1	3.00	2.00
Median	6.00	4.00
Q3	9.00	6.00
Max	47.00	38.00

4.3.3 Other categories

Chronic Disease Management, Nurse (10997), Medication Review (DMMR), Case Conference and Health Assessment (HA) billings
Period 1 - 12 month period prior to 2nd December 2019
Period 2 - 9 month period after and including 2nd December 2019

descr(active_patient_delta %>%
         select(CDM1, CDM2, Nurse1, Nurse2),
      stats = c("mean", "sd", "min", "q1", "med", "q3", "max"),
      style = "rmarkdown", headings = FALSE
)

	CDM1	CDM2	Nurse1	Nurse2
Mean	1.07	0.57	0.44	0.50
Std.Dev	1.34	0.97	0.89	0.97
Min	0.00	0.00	0.00	0.00
Q1	0.00	0.00	0.00	0.00
Median	0.00	0.00	0.00	0.00
Q3	2.00	1.00	1.00	1.00
Max	6.00	4.00	6.00	6.00

descr(x = (active_patient_delta %>%
             select(DMMR1, DMMR2, Conference1, Conference2, HA1, HA2)),
      stats = c("mean", "sd", "min", "q1", "med", "q3", "max"),
      style = "rmarkdown", headings = FALSE
)

	Conference1	Conference2	DMMR2	HA1	HA2
Mean	0.03	0.02	0.01	0.09	0.03
Std.Dev	0.20	0.15	0.12	0.29	0.16
Min	0.00	0.00	0.00	0.00	0.00
Q1	0.00	0.00	0.00	0.00	0.00
Median	0.00	0.00	0.00	0.00	0.00
Q3	0.00	0.00	0.00	0.00	0.00
Max	3.00	2.00	1.00	2.00	1.00

5 Estimations

The overall approach is a ‘difference-in-difference’ estimation.

compare the medication count and billings of patients in the time period when the pharmacist was working at the Kensington site (time period 2), with the same patients who were also seen at the clinic prior to the pharmacist commencing (time period 1).
Sex is also a (weak) predictor, and is included in the baseline model.
Some, but not all, of these same patients were also were contacted by the pharmacist, either in person or by telephone.

5.0.1 Analysis method

The choice of baseline model is described in Appendix A.

Subsequent analysis will be aimed at determining whether contact by the practice pharmacist influenced medications counts and billings in patients. The analysis will be determine whether contact by the pharmacist explains any of the ‘variance’ not explained in the ‘baseline model’ described in Appendix A.

5.1 Medication count

Contact by the pharmacist with a patient is strongly related to reduction in medication count (95% confidence interval -0.82 to -0.21, $p=0.001$).

MedCountGLM_base <- lm(
  MedCountDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
MedCountGLM <- lm(
  MedCountDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(MedCountGLM_base, MedCountGLM)

	MedCountDelta			MedCountDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	0.16	-0.03 – 0.36	0.096	0.14	-0.05 – 0.33	0.157
clinicVisit1	-0.04	-0.06 – -0.02	0.001	-0.03	-0.06 – -0.01	0.004
clinicVisit2	0.05	0.01 – 0.09	0.007	0.06	0.02 – 0.10	0.002
Sex [M]	0.12	-0.08 – 0.33	0.237	0.14	-0.06 – 0.35	0.171
DirectContactTRUE				-0.52	-0.83 – -0.22	0.001
Observations	1193			1193
R² / R² adjusted	0.011 / 0.009			0.021 / 0.017

5.2 Hospital admissions and discharge

The number of times a patient is seen in the clinic in time period 1 or 2, the number of medications and the number of admissions in time period 1 all are positively correlated with whether the patient was contacted/seen by the pharmacist in time period 2.

This suggests that the pharmacist was seeing patients who are at high risk, including at high risk of (re-)admission to hospital.

seen_by_pharmacist_model <- lm(
  pharmacistContact ~ clinicVisit1 + clinicVisit2 + MedCount1 + discharge1,
  data = active_patient_discharge
)
tab_model(seen_by_pharmacist_model, digits = 4)

	pharmacistContact
Predictors	Estimates	CI	p
(Intercept)	-0.0188	-0.0201 – -0.0175	<0.001
clinicVisit1	0.0016	0.0014 – 0.0017	<0.001
clinicVisit2	0.0039	0.0036 – 0.0042	<0.001
MedCount1	0.0038	0.0036 – 0.0040	<0.001
discharge1	0.0031	0.0023 – 0.0038	<0.001
Observations	137908
R² / R² adjusted	0.048 / 0.048

Even when comparing with ‘synthetic controls’, the patients seen by the pharmacist appear to be at high risk of being admitted to hospital.

The ‘synthetic’ controls are matched to patients seen by the pharmacist by the number of clinic visits in time period 1 and 2, the number of medications in the medication list at the end of time period 1, sex and age group (twenty year strata).

In the plots below, ‘week 53’ (indicated by a red dashed line) shows when the pharmacist contacted the patient. On visual inspection, there appears to be in a slight rise in the number of (pharmacist contact) ‘treatment’ patients being admitted to hospital when compared to the synthetic controls (right plot), even prior to week 53 (time of contact with the pharmacist).

df <- active_patient_discharge %>%
  filter(!is.na(Sex)) %>%
  mutate(Sex = dplyr::if_else(Sex == "F", 0, 1)) %>%
  as.data.frame(stringsasFactor = FALSE)

cov.var <- c("clinicVisit1", "clinicVisit2", "MedCount1", "Sex", "AgeGroup")

discharge_synth <- microsynth(data = df,
                              idvar = "InternalID", timevar = "timeVar",
                              intvar = "pharmacistContact",
                              match.out = list("admission" = c(2, 8, 8, 8)),
                              match.covar.min = cov.var,
                              result.var = "admission", omnibus.var = "admission",
                              n.cores = 1)
plot_microsynth(discharge_synth, main.diff = "admission trend vs control", xlab.tc = "weeks", xlab.diff = "weeks")

5.3 Billings

5.3.1 Total Billings and contact with pharmacist

Contact by the pharmacist with a patient is strongly related to additional Medicare billings in the selected categories (95% confidence interval 0.70 to 1.58 additional billings per patient, $p<0.001$).

TotalGLM_base <- lm(
  TotalDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
TotalGLM <- lm(
  TotalDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(TotalGLM_base, TotalGLM)

	TotalDelta			TotalDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	0.34	0.07 – 0.62	0.014	0.40	0.13 – 0.67	0.004
clinicVisit1	-0.95	-0.99 – -0.92	<0.001	-0.97	-1.00 – -0.94	<0.001
clinicVisit2	1.22	1.16 – 1.27	<0.001	1.20	1.14 – 1.25	<0.001
Sex [M]	-0.25	-0.54 – 0.04	0.093	-0.29	-0.58 – -0.00	0.047
DirectContactTRUE				1.14	0.70 – 1.58	<0.001
Observations	1193			1193
R² / R² adjusted	0.737 / 0.736			0.743 / 0.742

5.3.2 ‘Standard’ billings and contact with the pharmacist

Contact by the pharmacist with a patient is strongly related to additional ‘standard’ Medicare item number billings (95% confidence interval 0.48 to 1.25 additional billings per patient, $p<0.001$)

StandardGLM_base <- lm(
  StandardDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
StandardGLM <- lm(
  StandardDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(StandardGLM_base, StandardGLM)

	StandardDelta			StandardDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	0.54	0.30 – 0.77	<0.001	0.58	0.34 – 0.81	<0.001
clinicVisit1	-0.82	-0.85 – -0.80	<0.001	-0.84	-0.86 – -0.81	<0.001
clinicVisit2	1.02	0.98 – 1.07	<0.001	1.01	0.96 – 1.05	<0.001
Sex [M]	-0.13	-0.38 – 0.12	0.318	-0.16	-0.41 – 0.09	0.210
DirectContactTRUE				0.86	0.48 – 1.25	<0.001
Observations	1193			1193
R² / R² adjusted	0.732 / 0.731			0.736 / 0.735

5.3.3 Nurse provision of services related to chronic disease management

Contact by the pharmacist with a patient is strongly related to additional nurse chronic disease management item number billings (95% confidence interval 0.12 to 0.46, $p=0.001$).

NurseGLM_base <- lm(
  NurseDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
NurseGLM <- lm(
  NurseDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(NurseGLM_base, NurseGLM)

	NurseDelta			NurseDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	-0.03	-0.14 – 0.08	0.576	-0.02	-0.12 – 0.09	0.760
clinicVisit1	-0.04	-0.05 – -0.03	<0.001	-0.05	-0.06 – -0.03	<0.001
clinicVisit2	0.11	0.09 – 0.13	<0.001	0.11	0.09 – 0.13	<0.001
Sex [M]	-0.06	-0.17 – 0.06	0.320	-0.07	-0.18 – 0.05	0.238
DirectContactTRUE				0.29	0.12 – 0.46	0.001
Observations	1193			1193
R² / R² adjusted	0.085 / 0.082			0.093 / 0.090

5.3.4 Case conference item numbers

Contact by the pharmacist with a patient is weakly related to additional case conference item number claims (95% confidence interval -0.01 to 0.08, $p=0.136$).

ConferenceGLM_base <- lm(
  ConferenceDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
ConferenceGLM <- lm(
  ConferenceDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(ConferenceGLM_base, ConferenceGLM)

	ConferenceDelta			ConferenceDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	-0.02	-0.05 – 0.00	0.099	-0.02	-0.05 – 0.01	0.125
clinicVisit1	0.00	-0.00 – 0.00	0.399	0.00	-0.00 – 0.00	0.559
clinicVisit2	-0.00	-0.01 – 0.00	0.455	-0.00	-0.01 – 0.00	0.331
Sex [M]	0.01	-0.02 – 0.04	0.366	0.01	-0.02 – 0.04	0.412
DirectContactTRUE				0.03	-0.01 – 0.08	0.136
Observations	1193			1193
R² / R² adjusted	0.001 / -0.001			0.003 / -0.000

5.3.5 Domiciliary ‘Home’ medication reviews

Contact by a pharmacist with a patient is strongly related to home medication review claims (95% confidence interval 0.07 to 0.11, $p<0.001$)

DMMRGLM_base <- lm(
  DMMRDelta ~ clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
DMMRGLM <- lm(
  DMMRDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta
)
tab_model(DMMRGLM_base, DMMRGLM)

	DMMRDelta			DMMRDelta
Predictors	Estimates	CI	p	Estimates	CI	p
(Intercept)	-0.01	-0.03 – -0.00	0.033	-0.01	-0.02 – 0.00	0.134
clinicVisit1	-0.00	-0.00 – 0.00	0.827	-0.00	-0.00 – 0.00	0.081
clinicVisit2	0.01	0.00 – 0.01	<0.001	0.01	0.00 – 0.01	<0.001
Sex [M]	0.01	-0.01 – 0.02	0.248	0.00	-0.01 – 0.02	0.495
DirectContactTRUE				0.09	0.07 – 0.11	<0.001
Observations	1193			1193
R² / R² adjusted	0.039 / 0.036			0.102 / 0.099

6 Summary and Conclusion

During a period (December 2019 to August 2020) of considerable societal and medical practice upheaval due to a coronavirus outbreak, a pharmacist in practice was able to reduce the overall medication count (-0.52 per patient) and increase the billings (‘standard’ +0.86, ‘nursing’ +0.29, ‘case conference’ +0.03 and ‘medication review’ +0.09) of patients the pharmacist contacted.

Additional revenue per unique patient contact is estimated at $90.14 (Appendix B). One hundred and seventy unique patients were contacted, resulting in a total estimated revenue increase of $15323.80.

Although a pharmacist in practice was expected to reduce the medication count and increase the number of medication reviews, case conferences and perhaps other chronic disease management activities, the increase in other primary care team activity (‘standard’ item numbers, and ‘nurse’ related chronic disease management items) was even more pronounced. This is perhaps because the pharmacist played a role inside a primary care team, and the pharmacist’s interaction with the patient resulted in subsequent intentional activity with other primary care team members.

Although there was significant benefit, both medical and financial, with each pharmacist interaction with a patient, there was an overall struggle to arrange such an interaction in the first place. Initial plans to create these interactions, ranging from asthma plans and device training, screening for non-urgent review of chronic disease and daily screening of patient appointment lists for patients with poly-pharmacy were severely impacted by a global pandemic of a respiratory illness which resulted in far fewer patients willing to come physically to the clinic, or stay in the clinic for longer periods of time.

Only fairly late in the trial were programs successfully implemented to recruit patients to see the pharmacist without requiring initial physical presence in the clinic e.g. review of incoming and filed discharge correspondence.

Success of future pharmacy in practice programs will rely not only on the great medical and financial value which each interaction brings, but also attention to arranging such interactions and improved follow-through on potential medication review and care planning opportunities.

Appendix

A Baseline model

The baseline model estimates the number of billings (in select categories) a patient will have in the second period (after 2nd December 2020) based on the billings in the first period (prior to 2nd December 2020). The model bases the prediction on the number of times the patient visits the practice in-person (i.e. not telehealth or telephone consult).

Visits to the practice are used as a predictor, as this is also expected to be correlated to contact with the in-practice pharmacist. Until fairly late in the trial, strategies for patients to have a consult with the pharmacist mainly depended on patients physically in the practice being introduced (physically) to the pharmacist.

A.0.1 Billing change, baseline model

TotalDelta is the difference between Total2 (the total number of billings in chosen categories in the second time period) and Total1 (the total number of billings in the first time period) : Total2 - Total1.

The total number of billings in the first time period Total1 is expected to be positively correlated with the number of physical visits in the first time period clinicVisit1, and so we expect TotalDelta to be negatively correlated with clinicVisit1. This is confirmed in plot A of Figure A.1 TotalDelta vs. clinicVisit1.

total_clinic1 <- ggplot(
  active_patient_delta,
  aes(x = clinicVisit1, y = TotalDelta)
) + geom_point() + geom_smooth(method=lm)
total_clinic2 <- ggplot(
  active_patient_delta,
  aes(x = clinicVisit2, y = TotalDelta)
) + geom_point() + geom_smooth(method=lm)
cowplot::plot_grid(total_clinic1, total_clinic2,
                   labels = c("A", "B"))

Figure A.1: Change in billings vs. visits to clinic - period 1 and 2

Similarly, the total number of billings in the second time period Total2 is expected to be positively correlated with the number of physical visits in the second time period clinicVisit2. So it is expected TotalDelta will be positively correlated with clinicVisit2. Plot B of Figure A.1 ‘TotalDelta vs. clinicVisit2’ shows the relationship is not as strong/clear as in the case of TotalDelta vs clinicVisit1. This is most likely due to the presence of ‘non-physical’ billed contacts such as tele-health/telephone billed items.

Adding Sex to the model is of weak benefit. Adding DOB (date of birth) to the model did not provide any additional benefit.

clinicVisit1_model <- lm(TotalDelta ~ clinicVisit1, data = active_patient_delta)
clinicVisit2_model <- lm(TotalDelta ~ clinicVisit2, data = active_patient_delta)
clinicVisitBoth_model <- lm(TotalDelta ~ clinicVisit1 + clinicVisit2, data = active_patient_delta)
clinicVisitBothSex_model <- lm(TotalDelta ~ clinicVisit1 + clinicVisit2 + Sex, data = active_patient_delta)
tab_model(clinicVisit1_model, clinicVisit2_model, clinicVisitBoth_model, clinicVisitBothSex_model)

	TotalDelta			TotalDelta			TotalDelta			TotalDelta
Predictors	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p
(Intercept)	1.18	0.78 – 1.58	<0.001	-3.38	-3.79 – -2.97	<0.001	0.24	-0.01 – 0.49	0.056	0.34	0.07 – 0.62	0.014
clinicVisit1	-0.48	-0.52 – -0.44	<0.001				-0.95	-0.99 – -0.92	<0.001	-0.95	-0.99 – -0.92	<0.001
clinicVisit2				0.22	0.14 – 0.30	<0.001	1.22	1.17 – 1.27	<0.001	1.22	1.16 – 1.27	<0.001
Sex [M]										-0.25	-0.54 – 0.04	0.093
Observations	1196			1196			1196			1193
R² / R² adjusted	0.305 / 0.304			0.023 / 0.022			0.737 / 0.737			0.737 / 0.736

A model predicting TotalDelta by combining clinicVisit1, clinicVisit2 and Sex appears to be a reasonable model. Residuals of this model are shown in Figure A.2

par(mfrow = c(2,2))
plot(clinicVisitBothSex_model)

Figure A.2: Residuals of ‘TotalDelta ~ clinicVisit1 + clinicVisit2 + Sex’ model

A.0.2 Medication count change, baseline model

A similar model is build for changes in medication count MedCountDelta, depending on the number of physical visits in the two time periods (clinicVisit1 and clinicVisit2).

clinicVisit1 and clinicVisit2 do help predict MedCountDelta, but the relationship is not as strongs as for billings.

medcount_clinic1 <- ggplot(
  active_patient_delta,
  aes(x = clinicVisit1, y = MedCountDelta)
) + geom_point() + geom_smooth(method=lm)
medcount_clinic2 <- ggplot(
  active_patient_delta,
  aes(x = clinicVisit2, y = MedCountDelta)
) + geom_point() + geom_smooth(method=lm)
cowplot::plot_grid(medcount_clinic1, medcount_clinic2,
                   labels = c("A", "B"))

Figure A.3: Change in medication count vs. visits to clinic - period 1 and 2

medcount_clinicVisit1_model <- lm(MedCountDelta ~ clinicVisit1,
                                  data = active_patient_delta)
medcount_clinicVisit2_model <- lm(MedCountDelta ~ clinicVisit2,
                                  data = active_patient_delta)
medcount_clinicVisitBoth_model <- lm(MedCountDelta ~ clinicVisit1 + clinicVisit2,
                                     data = active_patient_delta)
medcount_clinicVisitBothSex_model <- lm(MedCountDelta ~ clinicVisit1 + clinicVisit2 + Sex,
                                        data = active_patient_delta)
tab_model(medcount_clinicVisit1_model, medcount_clinicVisit2_model,
          medcount_clinicVisitBoth_model, medcount_clinicVisitBothSex_model)

	MedCountDelta			MedCountDelta			MedCountDelta			MedCountDelta
Predictors	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p
(Intercept)	0.25	0.08 – 0.42	0.004	0.06	-0.09 – 0.21	0.447	0.21	0.04 – 0.39	0.016	0.16	-0.03 – 0.36	0.096
clinicVisit1	-0.02	-0.04 – -0.00	0.024				-0.04	-0.06 – -0.02	0.001	-0.04	-0.06 – -0.02	0.001
clinicVisit2				0.01	-0.02 – 0.04	0.554	0.05	0.01 – 0.09	0.008	0.05	0.01 – 0.09	0.007
Sex [M]										0.12	-0.08 – 0.33	0.237
Observations	1193			1193			1193			1193
R² / R² adjusted	0.004 / 0.003			0.000 / -0.001			0.010 / 0.009			0.011 / 0.009

par(mfrow = c(2,2))
plot(medcount_clinicVisitBothSex_model)

Figure A.4: Residuals of ‘MedCountDelta ~ clinicVisit1 + clinicVisit2 + Sex’ model

B Detailed revenue estimates

active_patient_billings_counts_specifics <- lapply(
  active_patient_billings_period,
  function(x) {
    # count the number of times each billing category occurred in each 9 month period
    x %>%
      group_by(InternalID) %>%
      summarise(StandardA = sum(MBSItem %in% c("3", "91790", "91795")),
                StandardB = sum(MBSItem %in% c("23", "91800", "91809")),
                StandardC = sum(MBSItem %in% c("36", "91801", "91810")),
                StandardD = sum(MBSItem %in% c("44", "91802", "91811")),
                StandardMH = sum(MBSItem %in% c("2713", "92115", "92127")),
                BBIncentive = sum(MBSItem %in% c("10990")),
                CDM721 = sum(MBSItem %in% c("721", "92024", "92068")),
                CDM723 = sum(MBSItem %in% c("723", "92025", "92069")),
                CDM732 = sum(MBSItem %in% c("732", "92028", "92072")),
                Conference735 = sum(MBSItem %in% c("735")),
                Conference739 = sum(MBSItem %in% c("739")),
                Conference743 = sum(MBSItem %in% c("743")),
                Conference747 = sum(MBSItem %in% c("747")),
                Conference750 = sum(MBSItem %in% c("750")),
                Conference758 = sum(MBSItem %in% c("758"))) %>%
      ungroup()           
    }
)

active_patient_delta_specifics <- active_patient_delta %>%
  left_join(active_patient_billings_counts_specifics[[1]] %>%
              rename(StandardA1 = StandardA,
                     StandardB1 = StandardB,
                     StandardC1 = StandardC,
                     StandardD1 = StandardD,
                     StandardMH1 = StandardMH,
                     BBIncentive1 = BBIncentive,
                     CDM721_1 = CDM721,
                     CDM723_1 = CDM723,
                     CDM732_1 = CDM732,
                     Conference735_1 = Conference735,
                     Conference739_1 = Conference739,
                     Conference743_1 = Conference743,
                     Conference747_1 = Conference747,
                     Conference750_1 = Conference750,
                     Conference758_1 = Conference758),
            by = "InternalID") %>%
  left_join(active_patient_billings_counts_specifics[[2]] %>%
              rename(StandardA2 = StandardA,
                     StandardB2 = StandardB,
                     StandardC2 = StandardC,
                     StandardD2 = StandardD,
                     StandardMH2 = StandardMH,
                     BBIncentive2 = BBIncentive,
                     CDM721_2 = CDM721,
                     CDM723_2 = CDM723,
                     CDM732_2 = CDM732,
                     Conference735_2 = Conference735,
                     Conference739_2 = Conference739,
                     Conference743_2 = Conference743,
                     Conference747_2 = Conference747,
                     Conference750_2 = Conference750,
                     Conference758_2 = Conference758),
            by = "InternalID") %>%
  replace_na(list(StandardA1 = 0, StandardA2 = 0,
                  StandardB1 = 0, StandardB2 = 0,
                  StandardC1 = 0, StandardC2 = 0,
                  StandardD1 = 0, StandardD2 = 0,
                  StandardMH1 = 0, StandardMH2 = 0,
                  BBIncentive1 = 0, BBIncentive2 = 0,
                  CDM721_1 = 0, CDM721_2 = 0,
                  CDM723_1 = 0, CDM723_2 = 0,
                  CDM732_1 = 0, CDM732_2 = 0,
                  Conference735_1 = 0, Confernece735_2 = 0,
                  Conference739_1 = 0, Confernece739_2 = 0,
                  Conference743_1 = 0, Confernece743_2 = 0,
                  Conference747_1 = 0, Confernece747_2 = 0,
                  Conference750_1 = 0, Confernece750_2 = 0,
                  Conference758_1 = 0, Confernece758_2 = 0)) %>%
  mutate(StandardADelta = StandardA2 - StandardA1,
         StandardBDelta = StandardB2 - StandardB1,
         StandardCDelta = StandardC2 - StandardC1,
         StandardDDelta = StandardD2 - StandardD1,
         BBIncentiveDelta = BBIncentive2 - BBIncentive1,
         CDM721Delta = CDM721_2 - CDM721_1,
         CDM723Delta = CDM723_2 - CDM723_1,
         CDM732Delta = CDM732_2 - CDM732_1,
         Conference735Delta = Conference735_2 - Conference735_1,
         Conference739Delta = Conference739_2 - Conference739_1,
         Conference743Delta = Conference743_2 - Conference743_1,
         Conference747Delta = Conference747_2 - Conference747_1,
         Conference750Delta = Conference750_2 - Conference750_1,
         Conference758Delta = Conference758_2 - Conference758_1)

B.1 ‘Standard’ item numbers

StandardAGLM <- lm(
  StandardADelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
StandardBGLM <- lm(
  StandardBDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
StandardCGLM <- lm(
  StandardCDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
StandardDGLM <- lm(
  StandardDDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
tab_model(StandardAGLM, StandardBGLM, StandardCGLM, StandardDGLM)

	StandardADelta			StandardBDelta			StandardCDelta			StandardDDelta
Predictors	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p
(Intercept)	0.14	-0.00 – 0.28	0.053	0.36	0.10 – 0.62	0.006	0.05	-0.15 – 0.26	0.601	-0.10	-0.23 – 0.03	0.149
DirectContactTRUE	0.01	-0.22 – 0.24	0.919	0.53	0.11 – 0.94	0.014	0.26	-0.07 – 0.59	0.121	0.10	-0.11 – 0.31	0.346
clinicVisit1	-0.08	-0.10 – -0.06	<0.001	-0.33	-0.36 – -0.30	<0.001	-0.31	-0.34 – -0.29	<0.001	-0.08	-0.09 – -0.06	<0.001
clinicVisit2	0.13	0.10 – 0.16	<0.001	0.41	0.36 – 0.47	<0.001	0.35	0.31 – 0.39	<0.001	0.09	0.07 – 0.12	<0.001
Sex [M]	-0.07	-0.22 – 0.08	0.370	-0.42	-0.69 – -0.14	0.003	0.21	-0.01 – 0.43	0.056	0.19	0.05 – 0.32	0.009
Observations	1193			1193			1193			1193
R² / R² adjusted	0.080 / 0.077			0.272 / 0.270			0.343 / 0.340			0.077 / 0.074

B.2 Bulk-bill incentives

BBIncentiveGLM <- lm(
  BBIncentiveDelta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
tab_model(BBIncentiveGLM)

	BBIncentiveDelta
Predictors	Estimates	CI	p
(Intercept)	0.39	0.02 – 0.76	0.037
DirectContactTRUE	1.28	0.69 – 1.87	<0.001
clinicVisit1	-0.81	-0.86 – -0.77	<0.001
clinicVisit2	0.91	0.84 – 0.99	<0.001
Sex [M]	-0.06	-0.44 – 0.33	0.779
Observations	1193
R² / R² adjusted	0.523 / 0.521

B.3 Case conference

Conference735GLM <- lm(
  Conference735Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
Conference739GLM <- lm(
  Conference739Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
Conference743GLM <- lm(
  Conference743Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
tab_model(Conference735GLM, Conference739GLM, Conference743GLM)

	Conference735Delta			Conference739Delta			Conference743Delta
Predictors	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p
(Intercept)	-0.00	-0.01 – 0.01	0.889	-0.00	-0.01 – 0.01	0.915	-0.00	-0.01 – 0.01	0.602
DirectContactTRUE	0.02	-0.01 – 0.04	0.146	0.03	0.01 – 0.05	0.002	0.00	-0.01 – 0.02	0.689
clinicVisit1	-0.00	-0.00 – -0.00	0.022	0.00	-0.00 – 0.00	0.346	0.00	-0.00 – 0.00	0.727
clinicVisit2	0.00	0.00 – 0.01	0.036	-0.00	-0.00 – 0.00	0.122	0.00	-0.00 – 0.00	0.547
Sex [M]	0.00	-0.01 – 0.02	0.533	-0.00	-0.01 – 0.01	0.686	0.00	-0.01 – 0.01	0.881
Observations	1193			1193			1193
R² / R² adjusted	0.007 / 0.004			0.010 / 0.007			0.002 / -0.002

Conference747GLM <- lm(
  Conference747Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
Conference750GLM <- lm(
  Conference750Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
Conference758GLM <- lm(
  Conference758Delta ~ DirectContact + clinicVisit1 + clinicVisit2 + Sex,
  data = active_patient_delta_specifics
)
tab_model(Conference747GLM, Conference750GLM, Conference758GLM)

	Conference747Delta			Conference750Delta			Conference758Delta
Predictors	Estimates	CI	p	Estimates	CI	p	Estimates	CI	p
(Intercept)	-0.01	-0.03 – 0.00	0.064	-0.00	-0.01 – 0.00	0.275	-0.00	-0.01 – 0.01	0.764
DirectContactTRUE	-0.00	-0.02 – 0.02	0.900	-0.01	-0.02 – 0.00	0.102	-0.00	-0.01 – 0.01	0.502
clinicVisit1	0.00	-0.00 – 0.00	0.241	0.00	-0.00 – 0.00	0.210	0.00	-0.00 – 0.00	0.225
clinicVisit2	-0.00	-0.01 – -0.00	0.026	0.00	-0.00 – 0.00	0.766	-0.00	-0.00 – -0.00	0.038
Sex [M]	0.01	-0.01 – 0.02	0.249	-0.00	-0.01 – 0.01	0.650	0.00	-0.00 – 0.01	0.401
Observations	1193			1193			1193
R² / R² adjusted	0.006 / 0.002			0.004 / 0.001			0.005 / 0.002

B.4 Totals

Estimated ‘Delta’ is per patient contacted.

tabulated_revenue <- data.frame(Item = character(), Delta = numeric(), ValuePerItem = numeric()) %>%
  add_row(Item = "A 'Item 3'", Delta = StandardAGLM$coefficients["DirectContactTRUE"], ValuePerItem = 17.75) %>%
  add_row(Item = "B 'Item 23'", Delta = StandardBGLM$coefficients["DirectContactTRUE"], ValuePerItem = 38.75) %>%
  add_row(Item = "C 'Item 36'", Delta = StandardCGLM$coefficients["DirectContactTRUE"], ValuePerItem = 75.05) %>%
  add_row(Item = "D 'Item 44'", Delta = StandardDGLM$coefficients["DirectContactTRUE"], ValuePerItem = 110.50) %>%
  add_row(Item = "BB Incentive '10990'", Delta = BBIncentiveGLM$coefficients["DirectContactTRUE"], ValuePerItem = 12.95) %>%
  add_row(Item = "Nurse CDM '10997'", Delta = NurseGLM$coefficients["DirectContactTRUE"], ValuePerItem = 12.40) %>%
  add_row(Item = "DMMR", Delta = DMMRGLM$coefficients["DirectContactTRUE"], ValuePerItem = 159.65) %>%
  add_row(Item = "Conference 735", Delta = Conference735GLM$coefficients["DirectContactTRUE"], ValuePerItem = 72.90) %>%
  add_row(Item = "Conference 739", Delta = Conference739GLM$coefficients["DirectContactTRUE"], ValuePerItem = 124.75) %>%
  add_row(Item = "Conference 743", Delta = Conference743GLM$coefficients["DirectContactTRUE"], ValuePerItem = 207.95) %>%
  add_row(Item = "Conference 747", Delta = Conference747GLM$coefficients["DirectContactTRUE"], ValuePerItem = 53.55) %>%
  add_row(Item = "Conference 750", Delta = Conference750GLM$coefficients["DirectContactTRUE"], ValuePerItem = 91.75) %>%
  add_row(Item = "Conference 758", Delta = Conference758GLM$coefficients["DirectContactTRUE"], ValuePerItem = 152.80) %>%
  mutate(TotalValue = Delta * ValuePerItem)
  
kable(tabulated_revenue, digits = 3) %>%
  kable_styling(bootstrap_options = c("striped"))

Item	Delta	ValuePerItem	TotalValue
A ‘Item 3’	0.012	17.75	0.212
B ‘Item 23’	0.525	38.75	20.346
C ‘Item 36’	0.262	75.05	19.670
D ‘Item 44’	0.100	110.50	11.090
BB Incentive ‘10990’	1.284	12.95	16.624
Nurse CDM ‘10997’	0.292	12.40	3.617
DMMR	0.092	159.65	14.745
Conference 735	0.016	72.90	1.134
Conference 739	0.029	124.75	3.649
Conference 743	0.003	207.95	0.529
Conference 747	-0.001	53.55	-0.074
Conference 750	-0.009	91.75	-0.871
Conference 758	-0.003	152.80	-0.532

Total additional revenue per patient contacted : $90.14

C Data extraction using dMeasure/GPstat

## Patients contacted by pharmacist
## uses DailyMeasure/GPstat! https://github.com/DavidPatShuiFong/DailyMeasure

library(pipeR)
library(dplyr)
library(dMeasure)

a <- dMeasure::dMeasure$new()
a$open_emr_db()
## active patients
kensington_clinicians <- a$UserConfig %>>%
  filter(
    grepl(
      "Kensington",
      purrr::map_chr(Location, function(x) paste(x, collapse = ", "))
      # map_chr will create a 'collapsed' version of all the
      # listed locations
    )
  )

active_patient <- a$list_contact_count(
  date_from = as.Date("2017-10-01"), date_to = as.Date("2019-09-30"),
  clinicians = kensington_clinicians$Fullname,
  min_contact = 3,
  contact_type = "Services",
  lazy = FALSE, store = FALSE
)

active_patient_ID <- active_patient %>>% dplyr::pull(InternalID)

active_patient_billings <- a$emr_db$conn() %>>%
  tbl("SERVICES") %>>%
  select("SERVICEID", "SERVICEDATE", "INVOICEID", "MBSITEM", "PAID", "DESCRIPTION") %>>%
  filter(INVOICEID %in% invoiceID,
         as.Date(ServiceDate) >= as.Date("2017-10-01")) %>>%
  collect() %>>%
  rename(InvoiceID = INVOICEID, ServiceDate = SERVICEDATE,
         MBSItem = MBSITEM, Paid = PAID, Description = DESCRIPTION ) %>>%
  mutate(Description = trimws(Description),
         ServiceDate = as.Date(ServiceDate)) %>>%
  left_join(invoices %>>% collect(),
            by = "InvoiceID")

#active_patient_billings <- a$db$services %>>%
#  filter(
#    InternalID %in% active_patient_ID,
#    ServiceDate >= as.Date("2017-10-01")
#  ) %>>% dplyr::collect()

PharmacistVisits <- a$db$visits %>% filter(DrName == "NameOfThePharmacist")
# InternalID, VisitType, VisitDate, UserID, DrName

## Process documents which include word 'discharge' in the Subject or Description
## Document '.csv' created with GPstat! https://github.com/DavidPatShuiFong/DailyMeasure

active_patient <- read.csv("active_patient.csv", stringsAsFactors = FALSE)

d <- read.csv("KensingtonDCsummary.csv", stringsAsFactors = FALSE) %>%
  select(InternalID, CorrespondenceDate, Category, Subject, Detail, Comment) %>%
  mutate(CorrespondenceDate = as.Date(CorrespondenceDate)) %>%
  filter(CorrespondenceDate >= as.Date("2019-12-02") - years(1)) %>%
  filter(InternalID %in% active_patient$InternalID) %>% # only active patients
  arrange(InternalID, CorrespondenceDate) %>%
  filter(!(grepl("clinic", Category, ignore.case = TRUE))) %>% # remove various letters
  filter(!(grepl("clinic", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("clinic", Comment, ignore.case = TRUE))) %>%
  filter(!(grepl("outpatient", Category, ignore.case = TRUE))) %>%
  filter(!(grepl("outpatient", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("outpatient", Comment, ignore.case = TRUE))) %>%
  filter(!(grepl("OPD", Category, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("OPD", Subject, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("OPD", Comment, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("FTA", Category, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("FTA", Subject, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("FTA", Comment, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("FTA", Detail, ignore.case = FALSE, fixed = TRUE))) %>%
  filter(!(grepl("did not attend", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("Specialist letter", Category, ignore.case = TRUE))) %>%
  filter(!(grepl("Specialist letter", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("Specialist letter", Comment, ignore.case = TRUE)))

d <- d %>%
  filter(!(grepl("encounter for food challenge", Detail, ignore.case = TRUE))) %>%
  filter(!(grepl("birth summary", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("obstetric discharge", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("livebirth", Comment, ignore.case = TRUE))) %>%
  filter(!(grepl("lscs", Comment, ignore.case = TRUE))) %>%
  filter(!(grepl("caesar", Comment, ignore.case = TRUE))) %>%
  filter(!(grepl("caeser", Comment, ignore.case = TRUE))) %>% # wrong spelling
  filter(!(grepl("caesar", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("transition", Subject, ignore.case = TRUE))) %>%
  filter(!(grepl("elective", Detail, ignore.case = TRUE))) %>%
  filter(!(grepl("elective", Comment, ignore.case = TRUE)))

  
d <- d %>%
  distinct(InternalID, CorrespondenceDate, .keep_all = TRUE) 
# remvve duplicate dates



d <- d %>%
  mutate(Proximate = FALSE) 
# is the entry 'proximate' to a recent entry?
# in which case it is likely to be multiple correspondence referring to the same episode

for (i in seq(nrow(d), 2)) {
  # work 'backwards' with 
  if (d[i, "InternalID"] == d[i-1, "InternalID"] &&
      d[i, "CorrespondenceDate"] - d[i - 1, "CorrespondenceDate"] < 5) {
    # same InternalID and less than 5 days difference
    d[i, "Proximate"] <- TRUE
  }
}

d <- d %>%
  filter(!Proximate) %>%
  select(-Proximate)