Skip to content

Instantly share code, notes, and snippets.

@sTeamTraen

sTeamTraen/Cadegiani-Table1s.R

Created October 2, 2021 19:20
Show Gist options
  • Save sTeamTraen/83d54a97c6a2ff4be2ef5a23c56bd7fc to your computer and use it in GitHub Desktop.
Save sTeamTraen/83d54a97c6a2ff4be2ef5a23c56bd7fc to your computer and use it in GitHub Desktop.
Download ZIP
Analysis of Table 1 from two (or three) articles by Cadegiani et al.
Raw
Cadegiani-Table1s.R
# Code to analyse Table 1 from some articles about the use of anti-androgens to treat Covid-19.
# By Nick Brown, October 2021.
# Licence: CC-0.
# See https://steamtraen.blogspot.com/2021/10/applying-john-carlisles-table-1.html
# Compute the chi-square p value for a given analysis modality (see the blog post).
mychi <- function (s1, s2, n1, n2, analysis)
{
mat <- matrix(c(s1, (n1 - s1), s2, (n2 - s2)), ncol=2)
no.zeroes <- substr(analysis, 2, 2) == "a"
if (no.zeroes && any(mat < 1)) {
return(NULL)
}
min.3 <- substr(analysis, 2, 2) == "c"
if (min.3 && any(mat < 3)) {
return(NULL)
}
analysis <- as.integer(substr(analysis, 1, 1))
use.fisher <- if (analysis == 1) 5 else if (analysis == 2) 0 else 9999999
func <- chisq.test
if (any(mat < use.fisher)) {
func <- fisher.test
}
p <- func(mat)$p.value
if (is.nan(p)) { # If two cells are zero, chisq.test() will return NaN as the p value (fisher.test gives 1.0).
p <- 1.0 # In that case, we force a p value of 1.0 (which we will later discard anyway by the 0.98 test).
}
return(p)
}
# Calculate the p value from a t test.
# we use rounding to make it the most favourable possible for the authors.
mytt <- function (m1, sd1, n1, m2, sd2, n2, dp=1)
{
if (dp > 0) {
sign <- if (m1 < m2) 1 else -1
delta <- (10 ^ (-dp -1)) * 5
m1 <- m1 + (delta * sign)
m2 <- m2 - (delta * sign)
sd1 <- sd1 + delta
sd2 <- sd2 + delta
}
se <- sqrt(
(
( ((n1 - 1) * (sd1 ^ 2))
+ ((n2 - 1) * (sd2 ^ 2))
)
/ (n1 + n2 - 2)
)
* ((1 / n1) + (1 / n2))
)
t <- (m1 - m2) / se
p <- 2 * (pt(abs(t), (n1 + n2 - 2), lower.tail=FALSE))
return (p)
}
# Set this to TRUE to see the tables of the first digits of the p values.
show.p.deciles <- FALSE
for (analysis in c("1a", "1b", "1c", "2a", "2b", "2c", "3a", "3b", "3c")) {
cat("Analysis ", analysis, "\n", sep="")
### Proxalutamide Significantly Accelerates Viral Clearance and Reduces Time to Clinical Remission in Patients with Mild to Moderate COVID-19 (10.7759/cureus.13492)
pv1 <- c()
#Gender
pv1 <- c(pv1, mychi(100, 71, 128, 108, analysis))
#Age (y/o)
pv1 <- c(pv1, mytt(44.5, 13.1, 171, 46.1, 12.7, 65))
#Time-to-treat
pv1 <- c(pv1, mytt(4.2, 1.8, 171, 4.0, 1.6, 65))
#Obesity
pv1 <- c(pv1, mychi(28, 9, 171, 65, analysis))
#Hypertension
pv1 <- c(pv1, mychi(39, 14, 171, 65, analysis))
#Myocardial infarction
pv1 <- c(pv1, mychi(2, 0, 171, 65, analysis))
#Stroke
pv1 <- c(pv1, mychi(2, 0, 171, 65, analysis))
#Congestive heart failure
pv1 <- c(pv1, mychi(1, 1, 171, 65, analysis))
#Lipid disorders
pv1 <- c(pv1, mychi(30, 16, 171, 65, analysis))
#DM (type 1 and 2)
pv1 <- c(pv1, mychi(15, 6, 171, 65, analysis))
#Prediabetes
pv1 <- c(pv1, mychi(18, 5, 171, 65, analysis))
#Dysglycemia
pv1 <- c(pv1, mychi(33, 11, 171, 65, analysis))
#Asthma
pv1 <- c(pv1, mychi(8, 4, 171, 65, analysis))
#COPD
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Chronic kidney disease
pv1 <- c(pv1, mychi(1, 0, 171, 65, analysis))
#Liver fibrosis/cirrhosis
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Clinical depression
pv1 <- c(pv1, mychi(11, 3, 171, 65, analysis))
#Anxiety
pv1 <- c(pv1, mychi(29, 7, 171, 65, analysis))
#Hypothyroidism
pv1 <- c(pv1, mychi(18, 7, 171, 65, analysis))
#Autoimmune disorders
pv1 <- c(pv1, mychi(1, 0, 171, 65, analysis))
#Cancer
pv1 <- c(pv1, mychi(3, 1, 171, 65, analysis))
#Beta-blocker
pv1 <- c(pv1, mychi(7, 3, 171, 65, analysis))
#ACEi
pv1 <- c(pv1, mychi(1, 0, 171, 65, analysis))
#ARB
pv1 <- c(pv1, mychi(25, 5, 171, 65, analysis))
#Loop diuretics
pv1 <- c(pv1, mychi(2, 0, 171, 65, analysis))
#Thiazide diuretics
pv1 <- c(pv1, mychi(7, 2, 171, 65, analysis))
#CCB
pv1 <- c(pv1, mychi(15, 3, 171, 65, analysis))
#Statins
pv1 <- c(pv1, mychi(13, 4, 171, 65, analysis))
#Others
pv1 <- c(pv1, mychi(2, 0, 171, 65, analysis))
#Aspirin
pv1 <- c(pv1, mychi(5, 1, 171, 65, analysis))
#Clopidogrel
pv1 <- c(pv1, mychi(0, 1 , 171, 65, analysis))
#Warfarin
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Xa factor inhibitors
pv1 <- c(pv1, mychi(1, 0, 171, 65, analysis))
#Direct thrombin inhibitors
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Heparin
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Metformin
pv1 <- c(pv1, mychi(16, 5, 171, 65, analysis))
#GLP1R analogs
pv1 <- c(pv1, mychi(10, 3, 171, 65, analysis))
#SGLT2 inhibitors
pv1 <- c(pv1, mychi(16, 3, 171, 65, analysis))
#DPP4 inhibitors
pv1 <- c(pv1, mychi(2, 1, 171, 65, analysis))
#Sulfonylureas
pv1 <- c(pv1, mychi(3, 0, 171, 65, analysis))
#Glitazones
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Acarbose
pv1 <- c(pv1, mychi(0, 0, 171, 65, analysis))
#Insulin
pv1 <- c(pv1, mychi(2, 0, 171, 65, analysis))
#Orlistat
pv1 <- c(pv1, mychi(3, 2, 171, 65, analysis))
#Levothyroxine
pv1 <- c(pv1, mychi(18, 7, 171, 65, analysis))
#Hypnotics
pv1 <- c(pv1, mychi(6, 1, 171, 65, analysis))
#SSRI
pv1 <- c(pv1, mychi(17, 4, 171, 65, analysis))
#Other antidepressants and mood stabilizers
pv1 <- c(pv1, mychi(8, 3, 171, 65, analysis))
#Benzodiazepines
pv1 <- c(pv1, mychi(4, 0, 171, 65, analysis))
#Atypical antipsychotics
pv1 <- c(pv1, mychi(8, 2, 171, 65, analysis))
#Omega-3
pv1 <- c(pv1, mychi(3, 0, 171, 65, analysis))
#Vitamin D
pv1 <- c(pv1, mychi(20, 5, 171, 65, analysis))
#Zinc
pv1 <- c(pv1, mychi(9, 2, 171, 65, analysis))
#Vitamin C
pv1 <- c(pv1, mychi(14, 5, 171, 65, analysis))
#Multivitamin
pv1 <- c(pv1, mychi(12, 4, 171, 65, analysis))
# Limit p values to a maximum of 0.98, to prevent "infinity" errors.
pv1[pv1 > 0.98] <- 0.98
# Now compute the overall probability by Stouffer's formula.
overall.pv1 <- 1 - pnorm(sum(sapply(pv1, qnorm))/sqrt(length(pv1)))
cat("Overall p value from \"Proxalutamide Significantly Accelerates Viral Clearance\" = ", overall.pv1, "\n", sep="")
if (show.p.deciles) {
cat("Table of first digits of p values\n")
cat(0:9, "\n", sep=" ")
cat(tabulate(floor(pv1 * 10) + 1), "\n", sep=" ")
}
cat("\n")
### Early Antiandrogen Therapy With Dutasteride Reduces Viral Shedding, Inflammatory Responses, and Time-to-Remission in Males With COVID-19 (10.7759/cureus.13047)
### Note: In the preprint (10.21203/rs.3.rs-135815/v1), the number of people in each condition was reversed (43 placebo, 44 treatment).
pv2 <- c()
#AGE (y/o)
pv2 <- c(pv2, mytt(43.8, 14.1, 44, 40, 10.8, 43))
#BMI (Kg/m2)
pv2 <- c(pv2, mytt(26.1, 2.2, 44, 26.0, 3.2, 43))
#Time-to-treat
pv2 <- c(pv2, mytt(4.4, 1.4, 44, 4.3, 1.4, 43))
#Obesity
pv2 <- c(pv2, mychi(4, 5, 44, 43, analysis))
#Hypertension
pv2 <- c(pv2, mychi(12, 7, 44, 43, analysis))
#Mi
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Stroke
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Hf
pv2 <- c(pv2, mychi(0, 1, 44, 43, analysis))
#Lipid disorders
pv2 <- c(pv2, mychi(9, 9, 44, 43, analysis))
#Other cardiac dysfunctions
pv2 <- c(pv2, mychi(0, 1, 44, 43, analysis))
#T1DM and T2DM
pv2 <- c(pv2, mychi(4, 4, 44, 43, analysis))
#Prediabetes
pv2 <- c(pv2, mychi(8, 5, 44, 43, analysis))
#Dysglycemia
pv2 <- c(pv2, mychi(12, 9, 44, 43, analysis))
#Asthma
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#COPD
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#CKD
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Liver fibrosis/cirrhosis
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#Clinical depression
pv2 <- c(pv2, mychi(2, 3, 44, 43, analysis))
#Anxiety
pv2 <- c(pv2, mychi(6, 8, 44, 43, analysis))
#ADHD
pv2 <- c(pv2, mychi(5, 3, 44, 43, analysis))
#Insomnia
pv2 <- c(pv2, mychi(4, 1, 44, 43, analysis))
#Hypothyroidism
pv2 <- c(pv2, mychi(1, 4, 44, 43, analysis))
#Autoimmune disorders
pv2 <- c(pv2, mychi(1, 1, 44, 43, analysis))
#Cancer
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Erectile dysfunction
pv2 <- c(pv2, mychi(2, 1, 44, 43, analysis))
#Hypogonadism
pv2 <- c(pv2, mychi(8, 8, 44, 43, analysis))
#BPH
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#Beta-blocker
pv2 <- c(pv2, mychi(2, 1, 44, 43, analysis))
#ECAi
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#ARB
pv2 <- c(pv2, mychi(11, 7, 44, 43, analysis))
#Loop diuretics
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#Thiazide diuretics
pv2 <- c(pv2, mychi(4, 1, 44, 43, analysis))
#CCB
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#Statins
pv2 <- c(pv2, mychi(8, 8, 44, 43, analysis))
#Others
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Aspirin
pv2 <- c(pv2, mychi(0, 1, 44, 43, analysis))
#Clopidogrel
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Warfarin
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Xa factor inhibitors
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Direct thrombin inhibitors
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Heparin
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Metformin
pv2 <- c(pv2, mychi(10, 6, 44, 43, analysis))
#GLP1R analogues
pv2 <- c(pv2, mychi(2, 3, 44, 43, analysis))
#SGLT2 inhibitors
pv2 <- c(pv2, mychi(7, 3, 44, 43, analysis))
#DPP4 inhibitors
pv2 <- c(pv2, mychi(2, 1, 44, 43, analysis))
#Sylfonylureas
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Glitazone
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Acarbose
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Insulin
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Orlistat
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#Levothyroxine
pv2 <- c(pv2, mychi(1, 4, 44, 43, analysis))
#Liothyronine
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#Testosterone
pv2 <- c(pv2, mychi(7, 6, 44, 43, analysis))
#Aromatase inhibitors or SERMs
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#Hypnotics
pv2 <- c(pv2, mychi(3, 1, 44, 43, analysis))
#Selective serotonin reuptaker inhibitors (SSRI)
pv2 <- c(pv2, mychi(3, 6, 44, 43, analysis))
#Other antidepressants and humor stabilizers
pv2 <- c(pv2, mychi(4, 4, 44, 43, analysis))
#Benzodiazepines
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#Atypical antipsychotics
pv2 <- c(pv2, mychi(3, 2, 44, 43, analysis))
#CNS stimulants
pv2 <- c(pv2, mychi(5, 5, 44, 43, analysis))
#Alpha-1 adrenaline blockers
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#GnRH analogues and inhibitors, NSAA and others antiandrogens
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Omega-3
pv2 <- c(pv2, mychi(3, 1, 44, 43, analysis))
#Vitamin D
pv2 <- c(pv2, mychi(13, 15, 44, 43, analysis))
#Zinc
pv2 <- c(pv2, mychi(6, 7, 44, 43, analysis))
#Biotin
pv2 <- c(pv2, mychi(1, 0, 44, 43, analysis))
#Vitamin C
pv2 <- c(pv2, mychi(8, 7, 44, 43, analysis))
#Multivitamin
pv2 <- c(pv2, mychi(1, 2, 44, 43, analysis))
#BCG vaccine
pv2 <- c(pv2, mychi(44, 43, 44, 43, analysis)) #corrected, as paper has too many people in the treatment group!
#Influeza vaccine (in 2020)
pv2 <- c(pv2, mychi(14, 14, 44, 43, analysis))
#Pneumococcal vaccine (since 2017)
pv2 <- c(pv2, mychi(3, 5, 44, 43, analysis))
#Current smoking
pv2 <- c(pv2, mychi(2, 2, 44, 43, analysis))
#Regular physical activity
pv2 <- c(pv2, mychi(28, 32, 44, 43, analysis))
#Ivermectin
pv2 <- c(pv2, mychi(6, 6, 44, 43, analysis))
#Hydroxychloroquine
pv2 <- c(pv2, mychi(4, 3, 44, 43, analysis))
#Xa factor inhibitors
pv2 <- c(pv2, mychi(8, 3, 44, 43, analysis))
#Enoxaparin
pv2 <- c(pv2, mychi(3, 3, 44, 43, analysis))
#Glucocorticoids
pv2 <- c(pv2, mychi(6, 7, 44, 43, analysis))
#Vitamin C
pv2 <- c(pv2, mychi(4, 4, 44, 43, analysis))
#Zinc
pv2 <- c(pv2, mychi(3, 3, 44, 43, analysis))
#Vitamin D
pv2 <- c(pv2, mychi(1, 2, 44, 43, analysis))
#Colchicine
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#Bromexhine
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
#N-acetylcysteine
pv2 <- c(pv2, mychi(0, 0, 44, 43, analysis))
# Limit p values to a maximum of 0.98, to prevent "infinity" errors.
pv2[pv2 > 0.98] <- 0.98
# Now compute the overall probability by Stouffer's formula.
overall.pv2 <- 1 - pnorm(sum(sapply(pv2, qnorm))/sqrt(length(pv2)))
cat("Overall p value from \"Early Antiandrogen Therapy With Dutasteride\" = ", overall.pv2, "\n", sep="")
if (show.p.deciles) {
cat("Table of first digits of p values\n")
cat(0:9, "\n", sep=" ")
cat(tabulate(floor(pv2 * 10) + 1), "\n", sep=" ")
}
cat("\n")
### 5-Alpha-Reductase Inhibitors Reduce Remission Time of COVID-19 (10.1101/2020.11.16.20232512)
pv3 <- c()
#age
pv3 <- c(pv3, mytt(42.0, 11.5, 66, 42.0, 12.6, 64))
#weight
pv3 <- c(pv3, mytt(85.5, 11.6, 66, 85.5, 9.8, 64))
#height
pv3 <- c(pv3, mytt(1.8, 0.1, 66, 1.8, 0.1, 64))
#BMI - not included here because it's redundant with height and weight.
# (If we use BMI instead of height/weight the result is to the authors' disadvantage.)
#pv3 <- c(pv3, mytt(26.5, 2.8, 66, 26.5, 3.2, 64))
#obesity
pv3 <- c(pv3, mychi(8, 10, 66, 64, analysis))
#Married?
pv3 <- c(pv3, mychi(47, 43, 66, 64, analysis))
#Live alone?
pv3 <- c(pv3, mychi(12, 13, 66, 64, analysis))
#Hypertension
pv3 <- c(pv3, mychi(18, 15, 66, 64, analysis))
#Mi
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Stroke
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Hf
pv3 <- c(pv3, mychi(0, 1, 66, 64, analysis))
#Lipid disorders
pv3 <- c(pv3, mychi(14, 16, 66, 64, analysis))
#Other cardiac dysfunctions
pv3 <- c(pv3, mychi(0, 3, 66, 64, analysis))
#Diabetes
pv3 <- c(pv3, mychi(7, 7, 66, 64, analysis))
#Pré-dm2
pv3 <- c(pv3, mychi(11, 9, 66, 64, analysis))
#Asthma
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#COPD
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Chronic renal dis.
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Liver fibrosis/chirrosis
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Clinical depression
pv3 <- c(pv3, mychi(5, 4, 66, 64, analysis))
#Anxiety
pv3 <- c(pv3, mychi(13, 12, 66, 64, analysis))
#ADHD
pv3 <- c(pv3, mychi(5, 6, 66, 64, analysis))
#Insomnia
pv3 <- c(pv3, mychi(4, 6, 66, 64, analysis))
#Hypothyroidism
pv3 <- c(pv3, mychi(3, 4, 66, 64, analysis))
#Autoimmune disorders
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Cancer
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Hypogonadism
pv3 <- c(pv3, mychi(13, 11, 66, 64, analysis))
#BPH
pv3 <- c(pv3, mychi(2, 3, 66, 64, analysis))
#Beta-blocker
pv3 <- c(pv3, mychi(3, 3, 66, 64, analysis))
#ECAi
pv3 <- c(pv3, mychi(4, 1, 66, 64, analysis))
#ARB
pv3 <- c(pv3, mychi(14, 14, 66, 64, analysis))
#Loop diur.
pv3 <- c(pv3, mychi(3, 0, 66, 64, analysis))
#Thiazide diuretics
pv3 <- c(pv3, mychi(4, 5, 66, 64, analysis))
#CCB
pv3 <- c(pv3, mychi(2, 4, 66, 64, analysis))
#Statins
pv3 <- c(pv3, mychi(14, 15, 66, 64, analysis))
#Others
pv3 <- c(pv3, mychi(0, 1, 66, 64, analysis))
#Aspirin
pv3 <- c(pv3, mychi(0, 2, 66, 64, analysis))
#Clopidogrel
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Warfarin
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Xa factor inhibitors
pv3 <- c(pv3, mychi(0, 1, 66, 64, analysis))
#Direct thrombin inhibitors
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Heparins
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Metformin
pv3 <- c(pv3, mychi(16, 14, 66, 64, analysis))
#GLP1R analogue
pv3 <- c(pv3, mychi(4, 8, 66, 64, analysis))
#SGLT2 inhibitors
pv3 <- c(pv3, mychi(8, 7, 66, 64, analysis))
#DPP4 inhibitors
pv3 <- c(pv3, mychi(3, 2, 66, 64, analysis))
#Sylfonylureas
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Glitazone
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Acarbose
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Insulin
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Orlistat
pv3 <- c(pv3, mychi(4, 5, 66, 64, analysis))
#Levothyroxine
pv3 <- c(pv3, mychi(3, 4, 66, 64, analysis))
#Liothyronine
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Testosterone
pv3 <- c(pv3, mychi(11, 9, 66, 64, analysis))
#Aromatase inhibitors or SERMs
pv3 <- c(pv3, mychi(2, 2, 66, 64, analysis))
#Hipnotics
pv3 <- c(pv3, mychi(4, 5, 66, 64, analysis))
#Selective serotonin reuptaker inhibitors (SSRI)
pv3 <- c(pv3, mychi(10, 9, 66, 64, analysis))
#Other antidepressants and humor stabilizers
pv3 <- c(pv3, mychi(7, 6, 66, 64, analysis))
#Benzodiazepines
pv3 <- c(pv3, mychi(2, 2, 66, 64, analysis))
#Atypical antipsychotics
pv3 <- c(pv3, mychi(4, 3, 66, 64, analysis))
#CNS stimulants
pv3 <- c(pv3, mychi(5, 8, 66, 64, analysis))
#Alpha-1 adren. Blockers
pv3 <- c(pv3, mychi(2, 3, 66, 64, analysis))
#Gnrh analogues and inh., NSAA, others
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
#Erectile dysfunction
pv3 <- c(pv3, mychi(2, 1, 66, 64, analysis))
#Omega-3
pv3 <- c(pv3, mychi(3, 1, 66, 64, analysis))
#Vitamin D
pv3 <- c(pv3, mychi(17, 17, 66, 64, analysis))
#Zinc
pv3 <- c(pv3, mychi(6, 7, 66, 64, analysis))
#Biotin
pv3 <- c(pv3, mychi(1, 0, 66, 64, analysis))
#Vitamin C
pv3 <- c(pv3, mychi(8, 8, 66, 64, analysis))
#Multivitamin
pv3 <- c(pv3, mychi(2, 2, 66, 64, analysis))
#BCG
pv3 <- c(pv3, mychi(64, 63, 66, 64, analysis))
#INFLUENZA (in 2020)
pv3 <- c(pv3, mychi(20, 21, 66, 64, analysis))
#PNEUMOCOCCAL (since 2017)
pv3 <- c(pv3, mychi(4, 8, 66, 64, analysis))
#CURRENT SMOKING
pv3 <- c(pv3, mychi(4, 6, 66, 64, analysis))
#Ivermectin
pv3 <- c(pv3, mychi(7, 6, 66, 64, analysis))
#Hydroxychlorouquine
pv3 <- c(pv3, mychi(4, 3, 66, 64, analysis))
#Xa FACTOR INHIBITORS
pv3 <- c(pv3, mychi(9, 6, 66, 64, analysis))
#Enoxaparin
pv3 <- c(pv3, mychi(4, 4, 66, 64, analysis))
#Glucocorticoids
pv3 <- c(pv3, mychi(7, 8, 66, 64, analysis))
#Vitamin c
pv3 <- c(pv3, mychi(4, 5, 66, 64, analysis))
#Zinc
pv3 <- c(pv3, mychi(3, 4, 66, 64, analysis))
#Vitamin d
pv3 <- c(pv3, mychi(1, 3, 66, 64, analysis))
#Colchicine
pv3 <- c(pv3, mychi(0, 0, 66, 64, analysis))
# Limit p values to a maximum of 0.98, to prevent "infinity" errors.
pv3[pv3 > 0.98] <- 0.98
# Now compute the overall probability by Stouffer's formula.
overall.pv3 <- 1 - pnorm(sum(sapply(pv3, qnorm))/sqrt(length(pv3)))
cat("Overall p value from \"5-Alpha-Reductase Inhibitors\" = ", overall.pv3, "\n", sep="")
if (show.p.deciles) {
cat("Table of first digits of p values\n")
cat(0:9, "\n", sep=" ")
cat(tabulate(floor(pv3 * 10) + 1), "\n", sep=" ")
}
cat("\n")
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment