# library(devtools) # dev_mode(on = T) # # install.packages("ggplot2") # withr::with_libpaths( # new = "C:/Rdevlibs", # install_github('tidyverse/ggplot2' # ,force = T)) # # sessionInfo() # library(ggplot2) # dev_mode(on = F) # update.packages(ask=F) # update.packages(checkBuilt = TRUE) #install.packages(c("curl", "httr")) # install.packages("mapdata") # install packages and load libraries # install.packages("jsonlite") # install.packages("countrycode") # install.packages("psych") # install.packages("beepr") # install.packages("ggmap") # install.packages("maptools") # install.packages("rworldmap") # install.packages("spam") # install.packages("devtools") # install.packages("compareGroups") # install.packages("XML") # ####install.packages("doParallel") # install.packages("doSNOW") # install.packages("dtplyr") # install.packages("lsr") # install.packages("gmodels") # install.packages("Rmisc") # source("http://www.rforge.net/NCStats/InstallNCStats.R") # install_github("kassambara/easyGgplot2") # install.packages("dplyr") # installXLSXsupport() # install.packages('lsmeans') # install.packages("svglite") # install.packages("scales") # install.packages("colorspace") # install.packages('haplo.stats') # install.packages("modEvA") # install.packages("pscl") # install.packages("AER") # install.packages("cowplot") # install.packages("survey") # install.packages("caret") # install.packages("ggplot2") # install.packages("lattice") # install.packages("maptools") # install.packages("stringi") # doInstall <- TRUE # toInstall <- c("maps", "ggplot2") # if(doInstall){install.packages(toInstall, repos = "http://cran.us.r-project.org")} # lapply(toInstall, library, character.only = TRUE) #install.packages("mice") #install.packages("ggalt") # install.packages("ggthemes") # library(ggthemes) # devtools::install_github("hrbrmstr/ggalt") # devtools::install_github("dkahle/ggmap") # install.packages("mapproj") #install.packages("classInt") library(classInt) #display.brewer.all() data(wrld_simpl) library("mgcv") library(lsmeans) library(Rmisc) library(gmodels) #install.packages("maps", repos = "http://cloud.r-project.org/") #install.packages("rworldmap", repos = "http://cloud.r-project.org/") library(rworldmap) library(countrycode) require(RColorBrewer) library(rgdal) library(lsr) library(dtplyr) library(devtools) library(ggplot2) require(reshape2) library(NCStats) library(zoo) library(plyr) library(reshape2) library(compareGroups) library(multcomp) library(easyGgplot2) library(beepr) library(jsonlite) library(countrycode) library(Hmisc) library(gridExtra) library(psych) library("ggmap") library(maptools) library(maps) library(rworldmap) library(dplyr) library(scales) library(MASS) library(gridExtra) library(grid) library(XML) library(modEvA) library(pscl) library(AER) # library(cowplot) library(survey) library(caret) library(lmtest) library(scales) library(beepr) library(ggalt) library(dismo) library(mapdata) library(mapproj) library(mice) library(ggplot2) library(maps) library(maptools) library(geosphere) geompointshape=21 geompointsize=1 statsize=8 xtextsize=18 xtitlesize=18 ytextsize=18 ytitlesize=18 plotlist<-list() give.n <- function(x){return(c(y = mean(x), label = length(x))) } give.mean.logscale <- function(x){return(c(y = mean(x), label = signif(10^(mean(x)),digits = 3) ))} give.sd.logscale <- function(x){return(c(y = mean(x), label = signif(10^(sd(x)),digits = 3) ))} give.meanse.logscale <- function(x){return(c(y = mean(x), label = paste( signif(10^(mean(x)),digits = 3), "+/-", signif(10^(1.96*se(x)),digits = 3), sep = "" ) ))} give.plusminus <- function(x){return(c(y = mean(x), label = "+/-" )) } give.n1 <- function(x){return(c(y = log10(2.5), label = length(x))) } give.mean1.logscale <- function(x){return(c(y = log10(40), label = signif(10^(mean(x)),digits = 3) ))} give.sd1.logscale <- function(x){return(c(y = log10(110), label = signif(10^(sd(x)),digits = 3) ))} # give.se.logscale <- function(x){return(c(y = mean(x), label = signif(10^(mean(x)),digits = 3))) } TheEpoch <- as.POSIXct("1970-01-01 00:00:00",format="%Y-%m-%d %H:%M:%S",tz="UTC") load("sugofinterest.RData") ######################### sugammadex:sugofinterest ################################## df = alljsonsplitbyuser[,c( "freq_final", "PROVIDERTYPE What is your level of medical training?", "maincountry_income_str", "APPIMPORTANCE Please rate the importance of the app to your practice:", "LENGTHOFPRACTICE How many years have you been in practice? (Not counting years in training)", "PRACTICEMODEL What is your anesthesia practice model?", "PRACTICETYPE What is your primary practice environment?", "PRACTICESIZE What are is your practice size?", "PRACTICECOMMUNITYTYPE What community does your practice primarily serve?", "SUGAMMADEX-02 Do you have access to sugammadex in your clinical practice?", "SUGAMMADEX-03 How many years have you been using sugammadex?", "SUGAMMADEX-04 How many times per week would estimate you use sugammadex?", "SUGAMMADEX-05 Which of the following statements is true?", "SUGAMMADEX-06 Which of the following best characterizes your use of sugammadex?", "SUGAMMADEX-07 Which of the following statements are true for your practice?", "SUGAMMADEX-08 Which of the following statements are true for your personal clinical practice?", "SUGAMMADEX-09 Which of the following adverse events, suspected or proven to be sugammadex related, have you observed in your personal clinical practice?", "SUGAMMADEX-10 Which of the following adverse events, suspected or proven to be sugammadex related, are you aware have occurred to a colleague in their clinical practice?" )] df$maincountry_income_str <- factor(df$maincountry_income_str, levels = c("Low income","Lower middle income", "Upper middle income","High income") , ordered = F ) colnames(df) <- c("Frequency", "Provider", "MainCountryIncome", "Importance", "LengthOfPractice", "PRACTICEMODEL", "PRACTICETYPE", "PRACTICESIZE", "PRACTICECOMMUNITYTYPE", "SUGAMMADEX02", "SUGAMMADEX03", "SUGAMMADEX04", "SUGAMMADEX05", "SUGAMMADEX06", "SUGAMMADEX07", "SUGAMMADEX08", "SUGAMMADEX09", "SUGAMMADEX10" ) ImportanceLevels = c("Absolutely Essential", "Very Important","Of Average Importance","Of Little Importance","Not Important At All") df$Importance <- factor(df$Importance, levels = ImportanceLevels , ordered = F) HCLevels = c( "Physician: Attending/Consultant", "Physician: Fellow/Resident/Registrar", "Anesthesia Assistant (PA)", "Nurse Anesthetist (CRNA)", "Anesthesia Technician", "Student AA", "Student Nurse Anesthetist", "Technically Trained in Anesthesia", "Respiratory Therapist" ) unique(df$Provider) df$Provider <- factor(df$Provider, levels = HCLevels # , # ordered = F ) ha <- list( Physician = c("Physician: Attending/Consultant"), PhysicianTrainee = c("Physician: Fellow/Resident/Registrar"), Anesthetist = c("Anesthesia Assistant (PA)","Nurse Anesthetist (CRNA)"), AnesthetistTrainee = c("Student AA","Student Nurse Anesthetist"), AnesthesiaTechnician = c("Anesthesia Technician"), TechnicallyTrainedinAnesthesia = c("Technically Trained in Anesthesia"), RespiratoryTherapist = c("Respiratory Therapist") ) for (i in 1:length(ha)) levels(df$Provider)[levels(df$Provider)%in%ha[[i]]] <- names(ha)[i] df$LengthOfPractice <- as.numeric(df$LengthOfPractice) # df$LengthOfPractice #bad data at 25 years (default) df$LengthOfPractice[which(df$LengthOfPractice == 25)] <- NA # PracLens = c("0-5 Years","6-10 Years","11-20 Years","21-30 Years",">30 Years") # df$PracticeLengthOfTime <- NA # df$PracticeLengthOfTime[which(df$LengthOfPractice <= 5)] <- PracLens[1] # df$PracticeLengthOfTime[which(df$LengthOfPractice > 5 & df$LengthOfPractice <= 10)] <- PracLens[2] # df$PracticeLengthOfTime[which(df$LengthOfPractice > 10 & df$LengthOfPractice <= 20)] <- PracLens[3] # df$PracticeLengthOfTime[which(df$LengthOfPractice > 20 & df$LengthOfPractice <= 30)] <- PracLens[4] # df$PracticeLengthOfTime[which(df$LengthOfPractice > 30)] <- PracLens[5] # df$PracticeLengthOfTime <- as.factor(df$PracticeLengthOfTime) # df$PracticeLengthOfTime<-(factor(df$PracticeLengthOfTime, levels = PracLens)) # levels(df$PracticeLengthOfTime) PracLens = c("0-5 Years","6-10 Years","11-20 Years","21-30 Years") df$PracticeLengthOfTime <- NA df$PracticeLengthOfTime[which(df$LengthOfPractice <= 5)] <- PracLens[1] df$PracticeLengthOfTime[which(df$LengthOfPractice > 5 & df$LengthOfPractice <= 10)] <- PracLens[2] df$PracticeLengthOfTime[which(df$LengthOfPractice > 10 & df$LengthOfPractice <= 20)] <- PracLens[3] df$PracticeLengthOfTime[which(df$LengthOfPractice > 20)] <- PracLens[4] df$PracticeLengthOfTime <- as.factor(df$PracticeLengthOfTime) df$PracticeLengthOfTime<-(factor(df$PracticeLengthOfTime, levels = PracLens)) # levels(df$PracticeLengthOfTime) df$PRACTICEMODEL<-as.factor(df$PRACTICEMODEL) #unique(df$PRACTICEMODEL) HCLevels = c( "Physician only", "Physician supervised, anesthesiologist on site", "Physician supervised, non-anesthesiologist physician on site", "Physician supervised, no physician on site", "No physician supervision", "Not an anesthesia provider" ) df$PRACTICEMODEL <- factor(df$PRACTICEMODEL, levels = HCLevels, ordered = F) # levels(df$PRACTICEMODEL) df$PRACTICETYPE<-as.factor(df$PRACTICETYPE) # unique(df$PRACTICETYPE) HCLevels = c( "Private clinic or office", "Local health clinic", "Ambulatory surgery center", "Small community hospital", "Large community hospital", "Academic department/University hospital" ) df$PRACTICETYPE <- factor(df$PRACTICETYPE, levels = HCLevels, ordered = F) # levels(df$PRACTICETYPE) df$PRACTICESIZE<-as.factor(df$PRACTICESIZE) # unique(df$PRACTICESIZE) HCLevels = c( "I am the only practitioner for large area", "One of several practitioners in the area", "Group practice 1-5 members", "Group 5-10 members", "Group 10-25 members", "Group 25-50 members", "Group > 50 members" ) df$PRACTICESIZE <- factor(df$PRACTICESIZE, levels = HCLevels, ordered = F) df$PracticeSizeGrouped <- df$PRACTICESIZE ha <- list( Solo = c( "I am the only practitioner for large area", "One of several practitioners in the area" ), SmallGroup_LT10 = c( "Group practice 1-5 members", "Group 5-10 members" ), MedGroup10_25 = c( "Group 10-25 members" ), LargeGroup_GT25 = c( "Group 25-50 members", "Group > 50 members" ) ) for (i in 1:length(ha)) levels(df$PracticeSizeGrouped)[levels(df$PracticeSizeGrouped)%in%ha[[i]]] <- names(ha)[i] table(df$PracticeSizeGrouped) df$PRACTICECOMMUNITYTYPE<-as.factor(df$PRACTICECOMMUNITYTYPE) # unique(df$PRACTICETYPE) HCLevels = c( "Private clinic or office", "Local health clinic", "Ambulatory surgery center", "Small community hospital", "Large community hospital", "Academic department/University hospital" ) df$PRACTICETYPE <- factor(df$PRACTICETYPE, levels = HCLevels, ordered = F) sugofinterest<-df sugofinterest$SUG03_ANSWERED <- 1 sugofinterest$SUG03_ANSWERED[which(is.na(sugofinterest$SUGAMMADEX03))] <- 0 sugofinterest$SUG10_ANSWERED <- 1 sugofinterest$SUG10_ANSWERED[which(is.na(sugofinterest$SUGAMMADEX10))] <- 0 sugofinterest$SUG0310_SUM <- sugofinterest$SUG10_ANSWERED + sugofinterest$SUG03_ANSWERED sugofinterest$SawFatigue <- NA sugofinterest$SawFatigue[which(sugofinterest$SUG0310_SUM == 1)] <- "Yes" sugofinterest$SawFatigue[which(sugofinterest$SUG0310_SUM == 2)] <- "No" sugofinterest$SawFatigue <- factor(sugofinterest$SawFatigue) t<-table(sugofinterest$SawFatigue) addmargins(t) prop.table(t) sugofinterest$LogFrequency <- log(sugofinterest$Frequency) summary(sugofinterest$LogFrequency) # sugofinterest$LogFrequencyAsFactor <-cut2(sugofinterest$LogFrequency,cuts) sugofinterest$LogFrequencyAsFactor <-cut(sugofinterest$LogFrequency,breaks=3) levels(sugofinterest$LogFrequencyAsFactor) table(df$SUGAMMADEX02, df$MainCountryIncome) signif(100*prop.table(table(df$SUGAMMADEX02, df$MainCountryIncome),margin = 2),digits=3) save(sugofinterest,file="sugofinterest.Rdata") ######################### sugammadex:SawFatigue vs Provider ################################## # p<-plot(SawFatigue~Frequency+Provider+MainCountryIncome+Importance+LengthOfPractice, # data = sugofinterest) fieldofint<-"Provider" m<-as.matrix(table(sugofinterest[,paste(fieldofint,sep = "")], sugofinterest$SawFatigue)) m2<-as.data.frame(m[,"Yes"]/(m[,"Yes"]+m[,"No"])) names(m2)<-c("SawFatigue") #names(m2)<-c("Percentage With Survey Fatigue") m2[,paste(fieldofint,sep = "")]<-rownames(m2) m3<-melt(m2,id.var=fieldofint) HCLevels = c( "Physician", "PhysicianTrainee", "Anesthetist", "AnesthetistTrainee", "AnesthesiaTechnician", "TechnicallyTrainedinAnesthesia", "RespiratoryTherapist" ) m3$Provider <- factor(m3$Provider, levels = HCLevels , ordered = TRUE ) ggplot(m3, aes(x=Provider, y=value, colour=variable, fill=Provider)) + geom_bar(stat="identity") + geom_text(color='black',aes(label = (m[,"Yes"]) ), vjust=2, size=10) + geom_text(color='black',aes(label = (m[,"Yes"]+m[,"No"]) ), vjust=4, size=10) + scale_y_continuous(labels=percent) + labs( y = "Fatigue Rate",x="" ) + theme(axis.text.x = element_text(angle=20, hjust=1, size=xtextsize), axis.text.y = element_text(size=ytextsize), axis.title.x = element_text(size=xtitlesize), axis.title.y = element_text(size=ytitlesize), legend.text = element_blank(), legend.position = "none", legend.title = element_blank() ) ggsave(paste("SawFatigue", fieldofint, ".jpeg")) addmargins(m) ######################### sugammadex:SawFatigue vs MainCountryIncome ################################## # p<-plot(SawFatigue~Frequency+Provider+MainCountryIncome+Importance+LengthOfPractice, # data = sugofinterest) fieldofint<-"MainCountryIncome" m<-as.matrix(table(sugofinterest[,paste(fieldofint,sep = "")], sugofinterest$SawFatigue)) m2<-as.data.frame(m[,"Yes"]/(m[,"Yes"]+m[,"No"])) names(m2)<-c("SawFatigue") #names(m2)<-c("Percentage With Survey Fatigue") m2[,paste(fieldofint,sep = "")]<-rownames(m2) m3<-melt(m2,id.var=fieldofint) m3$MainCountryIncome <- factor(m3$MainCountryIncome, levels = c("Low income","Lower middle income", "Upper middle income","High income"), ordered = TRUE) ggplot(m3, aes(x=MainCountryIncome, y=value, colour=variable, fill=MainCountryIncome)) + geom_bar(stat="identity") + geom_text(color='black',aes(label = (m[,"Yes"]) ), vjust=2, size=10) + geom_text(color='black',aes(label = (m[,"Yes"]+m[,"No"]) ), vjust=4, size=10) + scale_y_continuous(labels=percent) + labs( y = "Fatigue Rate", x="User's Primary Country:\n Designated World Bank Income Level" ) + theme(axis.text.x = element_text(angle=20, hjust=1, size=xtextsize), axis.text.y = element_text(size=ytextsize), axis.title.x = element_text(size=xtitlesize), axis.title.y = element_text(size=ytitlesize), legend.text = element_blank(), legend.position = "none", legend.title = element_blank() ) ggsave(paste("SawFatigue", fieldofint, ".jpeg")) addmargins(m) ######################### sugammadex:SawFatigue vs Importance ################################## # p<-plot(SawFatigue~Frequency+Provider+MainCountryIncome+Importance+LengthOfPractice, # data = sugofinterest) fieldofint<-"Importance" m<-as.matrix(table(sugofinterest[,paste(fieldofint,sep = "")], sugofinterest$SawFatigue)) m2<-as.data.frame(m[,"Yes"]/(m[,"Yes"]+m[,"No"])) names(m2)<-c("SawFatigue") #names(m2)<-c("Percentage With Survey Fatigue") m2[,paste(fieldofint,sep = "")]<-rownames(m2) m3<-melt(m2,id.var=fieldofint) ImportanceLevels = c("Absolutely Essential", "Very Important","Of Average Importance","Of Little Importance","Not Important At All") m3$Importance <- factor(m3$Importance, levels = ImportanceLevels, ordered = TRUE) ggplot(m3, aes(x=Importance, y=value, colour=variable, fill=Importance)) + geom_bar(stat="identity") + geom_text(color='black',aes(label = (m[,"Yes"]) ), vjust=2, size=10) + geom_text(color='black',aes(label = (m[,"Yes"]+m[,"No"]) ), vjust=4, size=10) + scale_y_continuous(labels=percent) + labs( y = "Fatigue Rate", x = "Importance of App to Personal Practice" ) + theme(axis.text.x = element_text(angle=20, hjust=0.75, size=xtextsize), axis.text.y = element_text(size=ytextsize), axis.title.x = element_text(size=xtitlesize), axis.title.y = element_text(size=ytitlesize), legend.text = element_blank(), legend.position = "none", legend.title = element_blank() ) ggsave(paste("SawFatigue", fieldofint, ".jpeg")) addmargins(m) ######################### sugammadex:SawFatigue vs LengthOfPractice ################################## # p<-plot(SawFatigue~Frequency+Provider+MainCountryIncome+Importance+LengthOfPractice, # data = sugofinterest) fieldofint<-"PracticeLengthOfTime" m<-as.matrix(table(sugofinterest[,paste(fieldofint,sep = "")], sugofinterest$SawFatigue)) m2<-as.data.frame(m[,"Yes"]/(m[,"Yes"]+m[,"No"])) names(m2)<-c("SawFatigue") #names(m2)<-c("Percentage With Survey Fatigue") m2[,paste(fieldofint,sep = "")]<-rownames(m2) m3<-melt(m2,id.var=fieldofint) # ImportanceLevels = c("Absolutely Essential", "Very Important","Of Average Importance","Of Little Importance","Not Important At All") # # m3$Importance <- factor(m3$Importance, # levels = ImportanceLevels, # ordered = TRUE) PracLens = c("0-5 Years","6-10 Years","11-20 Years","21-30 Years") m3$PracticeLengthOfTime<-(factor(m3$PracticeLengthOfTime, levels = PracLens)) ggplot(m3, aes(x=PracticeLengthOfTime, y=value, colour=variable, fill=PracticeLengthOfTime)) + geom_bar(stat="identity") + geom_text(color='black',aes(label = (m[,"Yes"]) ), vjust=2, size=10) + geom_text(color='black',aes(label = (m[,"Yes"]+m[,"No"]) ), vjust=4, size=10) + scale_y_continuous(labels=percent) + labs( y = "Fatigue Rate", x = "Length of Time in Practice" ) + theme(axis.text.x = element_text(angle=20, hjust=0.75, size=xtextsize), axis.text.y = element_text(size=ytextsize), axis.title.x = element_text(size=xtitlesize), axis.title.y = element_text(size=ytitlesize), legend.text = element_blank(), legend.position = "none", legend.title = element_blank() ) ggsave(paste("SawFatigue", fieldofint, ".jpeg")) addmargins(m) ######################### sugammadex:SawFatigue vs Frequency ################################## # p<-plot(SawFatigue~+Provider+MainCountryIncome+Importance+LengthOfPractice, # data = sugofinterest) fieldofint<-"LogFrequencyAsFactor" plot(SawFatigue~LogFrequencyAsFactor, data = sugofinterest) m<-as.matrix(table(sugofinterest[,paste(fieldofint,sep = "")], sugofinterest$SawFatigue)) m2<-as.data.frame(m[,"Yes"]/(m[,"Yes"]+m[,"No"])) names(m2)<-c("SawFatigue") m2[,paste(fieldofint,sep = "")]<-rownames(m2) m3<-melt(m2,id.var=fieldofint) ImportanceLevels = levels(sugofinterest$LogFrequencyAsFactor) m3$Importance <- factor(m3$Importance, levels = ImportanceLevels, ordered = TRUE) ggplot(m3, aes(x=LogFrequencyAsFactor, y=value, colour=variable, fill=LogFrequencyAsFactor)) + geom_point(stat="identity") + scale_y_continuous(labels=percent) + labs( y = "Fatigue Rate", x = "Normalized Frequency of App Use (Low -- High)" ) + theme(axis.text.x = element_blank(), axis.text.y = element_text(size=ytextsize), axis.title.x = element_text(size=xtitlesize), axis.title.y = element_text(size=ytitlesize), legend.text = element_blank(), legend.position = "none", legend.title = element_blank() ) + guides(fill=FALSE) ggsave(paste("SawFatigue", fieldofint, ".jpeg")) m3$seq<- seq(1:nrow (m3)) tlm<-lm (value~seq, data = m3) summary(tlm) addmargins(m) ######################### sugammadex:glm Provider ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ Provider, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$Provider ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$Provider [which(!is.na(sugofinterest$Provider ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"Provider") ######################### sugammadex:glm MainCountryIncome ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ MainCountryIncome, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$MainCountryIncome ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$MainCountryIncome [which(!is.na(sugofinterest$MainCountryIncome ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"MainCountryIncome") ######################### sugammadex:glm Importance ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ Importance, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$Importance ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$Importance [which(!is.na(sugofinterest$Importance ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"Importance") ######################### sugammadex:glm PracticeLengthOfTime ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ PracticeLengthOfTime, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$PracticeLengthOfTime ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$PracticeLengthOfTime [which(!is.na(sugofinterest$PracticeLengthOfTime ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"PracticeLengthOfTime") ######################### sugammadex:glm PRACTICEMODEL ################################## sugofinterest$PRACTICEMODEL # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ PRACTICEMODEL, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$PRACTICEMODEL ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$PRACTICEMODEL [which(!is.na(sugofinterest$PRACTICEMODEL ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"PRACTICEMODEL") ######################### sugammadex:glm PRACTICETYPE ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ PRACTICETYPE, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$PRACTICETYPE ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$PRACTICETYPE [which(!is.na(sugofinterest$PRACTICETYPE ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"PRACTICETYPE") ######################### sugammadex:glm MainCountryIncome ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ MainCountryIncome, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$MainCountryIncome ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$MainCountryIncome [which(!is.na(sugofinterest$MainCountryIncome ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"MainCountryIncome") ######################### sugammadex:glm PracticeSizeGrouped ################################## # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ PracticeSizeGrouped, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$PracticeSizeGrouped ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$PracticeSizeGrouped [which(!is.na(sugofinterest$PracticeSizeGrouped ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"PracticeSizeGrouped") ######################### sugammadex:glm PRACTICECOMMUNITYTYPE ################################## df$PRACTICECOMMUNITYTYPE # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ PRACTICECOMMUNITYTYPE, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$PRACTICECOMMUNITYTYPE ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$PRACTICECOMMUNITYTYPE [which(!is.na(sugofinterest$PRACTICECOMMUNITYTYPE ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"PRACTICECOMMUNITYTYPE") ######################### sugammadex:glm LogFrequencyAsFactor ################################## sugofinterest$LogFrequencyAsFactor # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ LogFrequencyAsFactor, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) # summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$LogFrequencyAsFactor ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$LogFrequencyAsFactor [which(!is.na(sugofinterest$LogFrequencyAsFactor ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"LogFrequencyAsFactor") ######################### sugammadex:glm ################################## qqnorm(sugofinterest$LogFrequency) # tglm<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance # +PracticeLengthOfTime, data = sugofinterest, family = "binomial") m1 <- glm( SawFatigue~ Frequency, family = "binomial", data=sugofinterest) m2 <- update(m1, . ~ . + 0) summary(m1) # summary(m2) #ests est <- as.data.frame(cbind(Estimate = coef(m2), confint(m2,level=0.95))) #sampsize (sampsize<-length(which(!is.na(sugofinterest$Frequency ) & !is.na(sugofinterest$SawFatigue) ))) #table tb<-as.data.frame( table(sugofinterest$Frequency [which(!is.na(sugofinterest$Frequency ) & !is.na(sugofinterest$SawFatigue) )])) estm1 <- as.data.frame(cbind(Estimate = coef(m1), confint(m1,level=0.95))) estm1 write.table(cbind(tb$Var1, tb$Freq, est, exp(estm1), 1/(1+1/(exp(est))), pValue=coef(summary(m1))[,4]),sep="\t") regTermTest(m1,"Frequency") ######################### sugammadex:glm Multiple ################################## m1<-glm(SawFatigue~Frequency+Provider+MainCountryIncome+Importance +PracticeLengthOfTime + PRACTICEMODEL + PRACTICETYPE + PRACTICECOMMUNITYTYPE + PracticeSizeGrouped, data = sugofinterest, family = "binomial") summary(m1) anova(m1,test="Chisq") regTermTest(m1,"Frequency") regTermTest(m1,"Provider") regTermTest(m1,"MainCountryIncome") regTermTest(m1,"Importance") regTermTest(m1,"PracticeLengthOfTime") regTermTest(m1,"PRACTICEMODEL") regTermTest(m1,"PRACTICETYPE") regTermTest(m1,"PRACTICECOMMUNITYTYPE") regTermTest(m1,"PracticeSizeGrouped") coef(tglm) table(sugofinterest$Provider,sugofinterest$SUG0310_ANSWERED) summary(sugofinterest$SUG03_ANSWERED) summary(sugofinterest$SUG10_ANSWERED) table( sugofinterest$Provider, sugofinterest$SUGAMMADEX02 ) suginc<-table( alljsonsplitbyuser$maincountry_income_str, alljsonsplitbyuser$`SUGAMMADEX-02 Do you have access to sugammadex in your clinical practice?`) csqtest<-chisq.test(suginc) chisqPostHoc(csqtest,digits = 10)