• PCR = polymerise chain reaction test for COVID-19
• Directly test for the presence of COVID-19 antigen (as opposed to testing for the bodies immune response)
• False negatives can occur up to 30% of the time
• Given as CT values (the number of cycles after which a PCR becomes positive, the lower the number, the more virus is present)
• NEG if \(\geq 35\)
• Measured at days 0-6 (D0, D1, D2, D3, D4, D5, D6)

raw.dat$D0

##  [1] 31  26  26  24  24  POS 28  POS POS POS POS ND  POS POS ND  POS 30  29  23 
## [20] 30  34  28  22  17  22  27  34  19  25  15  28  23  30  27  24  29 
## Levels: 15 17 19 22 23 24 25 26 27 28 29 30 31 34 ND POS

• Remove outcomes!

names(cleaned.dat)

##  [1] "Patient"             "Age"                 "Sex"                
##  [4] "ClinicalStatus"      "TimeBeforeInclusion" "Hydroxychloroquine" 
##  [7] "Concentration"       "ConcentrationTime"   "Azithromycin"       
## [10] "D0"                  "D1"                  "D2"                 
## [13] "D3"                  "D4"                  "D5"                 
## [16] "D6"

dat = cleaned.dat[,-c(11:16)]
names(dat)

##  [1] "Patient"             "Age"                 "Sex"                
##  [4] "ClinicalStatus"      "TimeBeforeInclusion" "Hydroxychloroquine" 
##  [7] "Concentration"       "ConcentrationTime"   "Azithromycin"       
## [10] "D0"

library("tidyverse")
dat = dat %>%
  mutate(Hydroxychloroquine = recode(Hydroxychloroquine, "Yes"= 1, "No"=0))

What do we want to know about treatment???

“French Confirmed COVID-19 patients were included in a single arm protocol from early March to March 16th, to receive 600mg of hydroxychloroquine daily and their viral load in nasopharyngeal swabs was tested daily in a hospital setting. Depending on their clinical presentation, azithromycin was added to the treatment. Untreated patients from another center and cases refusing the protocol were included as negative controls.”

“Patients who were proposed a treatment with hydroxychloroquine were recruited and managed in Marseille centre. Controls without hydroxychloroquine treatment were recruited in Marseille, Nice, Avignon and Briançon centers, all located in South France.”

=> Give up?

names(dat)

##  [1] "Patient"             "Age"                 "Sex"                
##  [4] "ClinicalStatus"      "TimeBeforeInclusion" "Hydroxychloroquine" 
##  [7] "Concentration"       "ConcentrationTime"   "Azithromycin"       
## [10] "D0"

• Age
• Sex
• ClinicalStatus
• TimeBeforeInclusion
• D0

Unconfounded???

boxplot(dat$Age~dat$Hydroxychloroquine)

table(dat$Hydroxychloroquine, dat$ClinicalStatus)

##    
##     Asymptomatic LRTI URTI
##   0            4    2   10
##   1            2    6   12

• URTI = Upper Respiratory Tract Infection
• LRTI = Lower Respiratory Tract Infection (Pneumonia confirmed in all)
• Again imbalanced, but again in favor of control

table(raw.dat$Hydroxychloroquine, !(raw.dat$D0 == "POS" | raw.dat$D0 == "ND"))

##      
##       FALSE TRUE
##   No     10    6
##   Yes     0   20

boxplot(as.numeric(as.character(raw.dat$D0))~raw.dat$Hydroxychloroquine)

## Warning in eval(predvars, data, env): NAs introduced by coercion

covs = c("Age", "Sex", "TimeBeforeInclusion", "Asymptomatic", "URTI","LRTI")
x = select(dat, one_of(covs))
cov.balance(x, dat$Hydroxychloroquine)

Error in plot.window(…) : need finite ‘xlim’ values

table(dat$Sex)

## 
##  F  M 
## 21 15

dat = mutate(dat, Female = recode(Sex, "F" = 1, "M" = 0))

“French Confirmed COVID-19 patients were included in a single arm protocol from early March to March 16th, to receive 600mg of hydroxychloroquine daily and their viral load in nasopharyngeal swabs was tested daily in a hospital setting. Depending on their clinical presentation, azithromycin was added to the treatment. Untreated patients from another center and cases refusing the protocol were included as negative controls.”

“Patients who were proposed a treatment with hydroxychloroquine were recruited and managed in Marseille centre. Controls without hydroxychloroquine treatment were recruited in Marseille, Nice, Avignon and Briançon centers, all located in South France.”

• TimeBeforeInclusion is based on onset of symptoms - missing for asymptomatic individuals.
  
• It’s well-balanced anyway, so exclude from propensity score model.)

ps.model = glm(Hydroxychloroquine~Age+Female+Asymptomatic+URTI, data = dat, family=binomial(link = "logit"))
summary(ps.model)$coefficients

##                 Estimate Std. Error    z value  Pr(>|z|)
## (Intercept)  -0.50966473  1.4482756 -0.3519114 0.7249047
## Age           0.03336289  0.0214164  1.5578195 0.1192760
## Female       -0.82264302  0.8115238 -1.0137017 0.3107251
## Asymptomatic -0.74480711  1.3941319 -0.5342444 0.5931725
## URTI         -0.23357741  1.0291615 -0.2269589 0.8204557

overlap(lps, dat$Hydroxychloroquine)

## [1] "6 controls below any treated"
## [1] "0 treated above any controls"

table(dat$Hydroxychloroquine)

## 
##  0  1 
## 16 20

• Should we trim???

dat[!kept,1:5]

##    Patient Age Sex ClinicalStatus TimeBeforeInclusion
## 1        1  10   M   Asymptomatic                  NA
## 2        2  12   F   Asymptomatic                  NA
## 3        3  14   F   Asymptomatic                  NA
## 4        4  10   M   Asymptomatic                  NA
## 12      12  23   F           URTI                   5
## 16      16  23   F           URTI                   6
## 18      18  54   M   Asymptomatic                  NA
## 20      20  59   F   Asymptomatic                  NA

• All the asymptomatic individuals (and 2 others)
• Easier to explain: just drop asymptomatic individuals?

• How many subclasses?

table(trimmed.dat$Hydroxychloroquine)

## 
##  0  1 
## 10 18

• How many subclasses?

table(dat$Hydroxychloroquine)

## 
##  0  1 
## 16 20

breaks = quantile(lps, c(.33, .67))
breaks

##        33%        67% 
## -0.0858544  0.5964168

subclass = subclasses(lps, breaks = breaks)
subclass

##  [1] 1 1 1 1 2 3 3 3 2 3 3 1 2 1 1 1 1 2 2 1 2 1 3 3 2 3 2 1 3 2 2 3 2 2 3 3

subclass.average(y, dat$Hydroxychloroquine, subclass = subclass)

## $`Difference in Means within Each Subclass`
## [1]  0.0000000 -0.6666667 -0.8333333
## 
## $`Weighted Average Difference in Means`
##      [,1]
## [1,] -0.5

• Exact same overall result.
• (Not surprising - balance on observed covariates was pretty good to begin with)
• With larger samples we might read into the differences across subclasses (treatment effect better for those most likely to be treated), but within subclasses samples too small to say much of anything here

sum(is.na(y))

## [1] 6

table(is.na(y), dat$Hydroxychloroquine)

##        
##          0  1
##   FALSE 11 19
##   TRUE   5  1

• Missing at random?
• How to deal with a variable that’s part numeric and part factor in imputation???

dat = raw.dat
for (j in 10:16) dat[,j][dat[,j]=="ND"] = NA
dat = droplevels(dat)
M = 5
imp = mice(dat, m = M)
imp.dat = vector("list", M)
for (m in 1:M) imp.dat[[m]] = complete(imp, action = m)
save(imp.dat, file="imp.rData")

M=5
load("imp.rData")

Imputing after dichotomizing PCR results:

dat = cleaned.dat
imp = mice(dat, m = M)
imp.dat = vector("list", M)
for (m in 1:M) imp.dat[[m]] = complete(imp, action = m)
for (m in 1:M) ests[m] = mean(imp.dat[[m]]$D6[imp.dat[[m]]$Hydroxychloroquine=="Yes"]=="POS")-mean(imp.dat[[m]]$D6[imp.dat[[m]]$Hydroxychloroquine=="No"]=="POS")

ests

## [1] -0.575 -0.450 -0.575 -0.450 -0.575

mean(ests)

## [1] -0.525

• We can only impute the missing outcomes for the people that remained in the study but had “ND” PCR tests.

• However, 6 other units attrited, whom we have no data for.

• Nothing we can do about them analytically

• But we should investigate reasons for attrition

“Six hydroxychloroquine-treated patients were lost in follow-up during the survey because of early cessation of treatment. Reasons are as follows: three patients were transferred to intensive care unit, including one transferred on day2 post-inclusion who was PCR-positive on day1, one transferred on day3 post-inclusion who was PCR-positive on days1-2 and one transferred on day4 post-inclusion who was PCR-positive on day1 and day3; one patient died on day3 post inclusion and was PCR-negative on day2; one patient decided to leave the hospital on day3 post-inclusion and was PCR-negative on days1-2; finally, one patient stopped the treatment on day3 post-inclusion because of nausea and was PCR-positive on days1-2-3. The results presented here are therefore those of 36 patients (20 hydroxychloroquine-treated patients and 16 control patients). None of the control patients was lost in follow-up.”

THIS IS A HUGE PROBLEM!!!

• Azithromycin is an antibiotic

• “Depending on their clinical presentation, azithromycin was added to the treatment.”

table(dat$Hydroxychloroquine, dat$Azithromycin)

##      
##       No Yes
##   No  16   0
##   Yes 14   6

• Is the effect due to hydroxychloroquine or azithromycin??

• We can look at the units who received only hydroxychloroquine, but this is injecting certain bias, because patients only received this as needed (probably the sickest patients)

prop.table(table(dat$Hydroxychloroquine[dat$Azithromycin=="No"], y[dat$Azithromycin=="No"]),1)

##      
##               0         1
##   No  0.1818182 0.8181818
##   Yes 0.5384615 0.4615385

• Estimated treatment effect: -0.36

Excluding those who received Azithromycin:

prop.test(table(dat$Hydroxychloroquine[dat$Azithromycin=="No"], y[dat$Azithromycin=="No"]))

## Warning in prop.test(table(dat$Hydroxychloroquine[dat$Azithromycin == "No"], :
## Chi-squared approximation may be incorrect

## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  table(dat$Hydroxychloroquine[dat$Azithromycin == "No"], y[dat$Azithromycin ==     "No"])
## X-squared = 1.8909, df = 1, p-value = 0.1691
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.79466068  0.08137397
## sample estimates:
##    prop 1    prop 2 
## 0.1818182 0.5384615

We are now excluding the sickest units from the treatment group for several reasons:

• got transferred to ICU
• died
• needed antibiotics

It’s not surprising that those left are healthier!

“Hospitalized patients with confirmed COVID-19 were included in this study if they fulfilled two primary criteria: i) age >12 years; ii) PCR documented SARS-CoV-2 carriage in nasopharyngeal sample at admission whatever their clinical status.”

min(dat$Age)

## [1] 10

dat$D0

##  [1] POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  <NA> POS  POS  <NA>
## [16] POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS  POS 
## [31] POS  POS  POS  POS  POS  POS 
## Levels: POS

“In the EU Clinical Trial Register page, the study was described as evaluating PCR data on Day 1, Day 4, Day 7 and Day 14. However, the study show the data for Day 6, which is different than planned. Why did the authors not show the results on Day 7? Did the data on day 7 not look as good?”

“Assuming a 50% efficacy of hydroxychloroquine in reducing the viral load at day 7, a 85% power, a type I error rate of 5% and 10% loss to follow-up, we calculated that a total of 48 COVID-19 patients (ie, 24 cases in the hydroxychloroquine group and 24 in the control group) would be required for the analysis (Fleiss with CC).”

table(dat$Hydroxychloroquine)

## 
##  No Yes 
##  16  20

Even with the 6 missing, that still only brings us up to 42…

Hydroxychloroquine Update For April 6

Summarizes several other studies…