Hi William,
Thanks, I'll give that a shot. I tried using withCallingHandlers without
success but II admit I'm not familiar with it and may have used it wrong.
I'll report back.
Jen
On Tue, Mar 6, 2018, 5:42 PM William Dunlap <wdun...@tibco.com> wrote:
> You can capture warnings by using withCallingHandlers. Here is an
> example,
> its help file has more information.
>
> dataList <- list(
> A = data.frame(y=c(TRUE,TRUE,TRUE,FALSE,FALSE), x=1:5),
> B = data.frame(y=c(TRUE,TRUE,FALSE,TRUE,FALSE), x=1:5),
> C = data.frame(y=c(FALSE,FALSE,TRUE,TRUE,TRUE), x=1:5))
>
> withWarnings <- function(expr) {
> .warnings <- NULL # warning handler will append to this using '<<-'
> value <- withCallingHandlers(expr,
> warning=function(e) {
> .warnings <<- c(.warnings,
> conditionMessage(e))
> invokeRestart("muffleWarning")
> })
> structure(value, warnings=.warnings)
> }
> z <- lapply(dataList, function(data) withWarnings(coef(glm(data=data, y ~
> x, family=binomial))))
> z
>
> The last line produces
>
> > z
> $A
> (Intercept) x
> 160.80782 -45.97184
> attr(,"warnings")
> [1] "glm.fit: fitted probabilities numerically 0 or 1 occurred"
>
> $B
> (Intercept) x
> 3.893967 -1.090426
>
> $C
> (Intercept) x
> -115.02321 45.97184
> attr(,"warnings")
> [1] "glm.fit: fitted probabilities numerically 0 or 1 occurred"
>
> and lapply(z, attr, "warnings") will give you the warnings themselves.
>
>
>
> Bill Dunlap
> TIBCO Software
> wdunlap tibco.com
>
> On Tue, Mar 6, 2018 at 2:26 PM, Jen <plessthanpointohf...@gmail.com>
> wrote:
>
>> Hi, I am trying to automate the creation of tables for some simply
>> analyses. There are lots and lots of tables, thus the creation of a
>> user-defined function to make and output them to excel.
>>
>> My problem is that some of the analyses have convergence issues, which I
>> want captured and included in the output so the folks looking at them know
>> how to view those estimates.
>>
>> I am successfully able to do this in a straightforward set of steps.
>> However, once I place those steps inside a function it fails.
>>
>> Here's the code (sorry this is a long post):
>>
>> # create data
>> wt <- rgamma(6065, 0.7057511981, 0.0005502062)
>> grp <- sample(c(replicate(315, "Group1"), replicate(3672, "Group2"),
>> replicate(1080, "Group3"), replicate(998, "Group4")))
>> dta <- data.frame(grp, wt)
>> head(dta)
>> str(dta)
>>
>> # declare design
>> my.svy <- svydesign(ids=~1, weights=~wt, data=dta)
>>
>> # subset
>> grp1 <- subset(my.svy, grp == "Group1")
>>
>> # set options and clear old warnings
>> options(warn=0)
>> assign("last.warning", NULL, envir = baseenv())
>>
>> ## proportions and CIs
>> p <- ((svyciprop(~grp, grp1, family=quasibinomial))[1])
>>
>> # save warnings
>> wrn1 <- warnings(p)
>>
>> ci_l <- (confint(svyciprop(~grp, grp1, family=quasibinomial), 'ci')[1])
>> ci_u <- (confint(svyciprop(~grp, grp1, family=quasibinomial), 'ci')[2])
>>
>> ## sample counts
>> n <- unwtd.count(~grp, grp1)[1]
>>
>> ## combine into table
>> overall <- data.frame(n, p, ci_l, ci_u)
>> colnames(overall) <- c("counts", "Group1", "LL", "UL")
>>
>> ## add any warnings
>> ind <- length(wrn1)
>> ind
>>
>> if (ind == 0) { msg <- "No warnings" }
>> if (ind > 0) {msg <- names(warnings()) }
>> overall[1,5] <- msg
>>
>> print(overall)
>>
>> Here's the output from the above:
>>
>> > # set options and clear old warnings
>> > options(warn=0)
>> > assign("last.warning", NULL, envir = baseenv())
>> >
>> > ## proportions and CIs
>> > p <- ((svyciprop(~grp, grp1, family=quasibinomial))[1])
>> Warning message:
>> glm.fit: algorithm did not converge
>> >
>> > # save warnings
>> > wrn1 <- warnings(p)
>> >
>> > ci_l <- (confint(svyciprop(~grp, grp1, family=quasibinomial), 'ci')[1])
>> Warning message:
>> glm.fit: algorithm did not converge
>> > ci_u <- (confint(svyciprop(~grp, grp1, family=quasibinomial), 'ci')[2])
>> Warning message:
>> glm.fit: algorithm did not converge
>> >
>> > ## sample counts
>> > n <- unwtd.count(~grp, grp1)[1]
>> >
>> > ## combine into table
>> > overall <- data.frame(n, p, ci_l, ci_u)
>> > colnames(overall) <- c("counts", "Group1", "LL", "UL")
>> >
>> > ## add any warnings
>> > ind <- length(wrn1)
>> > ind
>> [1] 1
>> >
>> > if (ind == 0) { msg <- "No warnings" }
>> > if (ind > 0) {msg <- names(warnings()) }
>> > overall[1,5] <- msg
>> >
>> > print(overall)
>> counts Group1 LL UL
>> V5
>> counts 315 2.364636e-12 2.002372e-12 2.792441e-12 glm.fit: algorithm
>> did
>> not converge
>>
>> Here's the function:
>>
>> est <- function(var) {
>>
>> ## set up formula
>> formula <- paste ("~", var)
>>
>> ## set options and clear old warning
>> options(warn=0)
>> assign("last.warning", NULL, envir = baseenv())
>>
>> ## proportions and CIs
>> p <- ((svyciprop(as.formula(formula), grp1, family=quasibinomial))[1])
>>
>> ## save warnings
>> wrn1 <- warnings(p)
>>
>> ci_l <- (confint(svyciprop(as.formula(formula) , grp1,
>> family=quasibinomial), 'ci')[1])
>> ci_u <- (confint(svyciprop(as.formula(formula) , grp1,
>> family=quasibinomial), 'ci')[2])
>>
>> ## sample counts
>> n <- unwtd.count(as.formula(formula), grp1)[1]
>>
>> ## combine into table
>> overall <- data.frame(n, p, ci_l, ci_u)
>> colnames(overall) <- c("counts", "Group1", "LL", "UL")
>>
>>
>> ## add any warnings
>> ind <- length(warnings(p))
>> print(ind)
>>
>> if (ind == 0) { msg <- "No warnings" }
>> if (ind > 0) {msg <- names(warnings()) }
>> overall[1,5] <- msg
>>
>> print(overall)
>>
>> }
>>
>> Here's the output from running the function:
>>
>> > est("grp")
>> [1] 0
>> counts Group1 LL UL V5
>> counts 315 2.364636e-12 2.002372e-12 2.792441e-12 No warnings
>> Warning messages:
>> 1: glm.fit: algorithm did not converge
>> 2: glm.fit: algorithm did not converge
>> 3: glm.fit: algorithm did not converge
>>
>> So, the warnings are showing up in the output at the end of the function
>> but they're not being captured like they are when run outside of the
>> function. Note the 0 output from print(ind) and V7 has "No warnings".
>> I know a lot of things "behave" differently inside functions. Case in
>> point, the use of "as.formula(var)" rather than just "~grp" being passed
>> to
>> the function.
>>
>> I've failed to find a solution after much searching of various R related
>> forums. I even posted this to stackoverflow but with no response. So, if
>> anyone can help, I'd be appreciative.
>>
>> (sidenote: I used rgamma to create my sampling weights because that's what
>> most resembles the distribution of my weights and it's close enough to
>> reproduce the convergence issue. If I used rnorm or even rlnorm or
>> rweibull
>> I couldn't reproduce it. Just FYI.)
>>
>> Best,
>>
>> Jen
>>
>> [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
