[R] Visualization of coefficients
Achim Zeileis
Achim.Zeileis at uibk.ac.at
Wed Jul 7 11:31:06 CEST 2010
On Wed, 7 Jul 2010, Tal Galili wrote:
> Hi Achim and Allan,I updated the post with Allan's example (thanks Allan).
Thanks!
> Achim, you wrote:
> "Finally, the Poisson model in comparison with the binomial models does not
> make much sense, I guess."
> I agree. I wanted something to showcase the function on 3 models (with the
> same predictors), and that's the easiest I could think of. If you'd think
> of a smarter example I'd be happy to incorporate it.
You could generate Poisson data and then fit the binomial model to the
threshold version of the response.
But I guess that would be a bit over the top. Also, one could argue in
that case that a complementary log-log link should be employed.
Hence, I would simply "say" (verbally) that it works for baysglm, glm, lm,
polr objects and that a default method is available which takes
pre-computed coefficients and associated standard errors from any suitable
model.
Best,
Z
> Best,
> Tal
>
>
>
> ----------------Contact
> Details:-------------------------------------------------------
> Contact me: Tal.Galili at gmail.com | 972-52-7275845
> Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
> www.r-statistics.com (English)
> ---------------------------------------------------------------------------
> -------------------
>
>
>
>
> On Wed, Jul 7, 2010 at 12:10 PM, Achim Zeileis <Achim.Zeileis at uibk.ac.at>
> wrote:
> On Wed, 7 Jul 2010, Tal Galili wrote:
>
> Hello David,
> Thanks to your posting I started looking at the
> function in the arm package.
> It appears this function is quite mature, and
> offers (for example) the
> ability to easily overlap coefficients from several
> models.
>
>
> Re: more mature. arm's coefplot() is more flexible in certain
> respects, mine is more convenient in others. The overlay functionality
> is something arm's coefplot() is better in and it also as some further
> options (vertical vs. horizontal etc.). My coefplot() has the
> advantage that it does not need any modification as long as coef() and
> vcov() methods are available. Furthermore, "level" can specify the
> significance level (instead of always using one and two standard
> errors, respectively).
> But it shouldn't be too hard to create a superset of all options.
>
> I updated the post I published on the subject, so at the
> end of it I give an
> example of comparing the coef of several models:
> http://www.r-statistics.com/2010/07/visualization-of-regression-coefficient
>
> s-in-r/
>
>
> As Allan pointed out in his reply, something fully reproducible would
> be nice. Also, if you keep the example with quasi-complete separation,
> it would be worth pointing this out. (Because the maximum likelihood
> estimator is Infinity in this case.)
>
> Finally, the Poisson model in comparison with the binomial models does
> not make much sense, I guess.
>
> Best,
> Z
>
> Thanks again for the pointer.
>
> Best,
> Tal
>
>
>
>
> ----------------Contact
> Details:-------------------------------------------------------
> Contact me: Tal.Galili at gmail.com | 972-52-7275845
> Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il
> (Hebrew) |
> www.r-statistics.com (English)
> ---------------------------------------------------------------------------
>
> -------------------
>
>
>
>
> On Wed, Jul 7, 2010 at 12:02 AM, David Atkins
> <datkins at u.washington.edu>
> wrote:
>
>
> FYI, there is already a function coefplot in the arm
> package;
> for example, compare:
>
> > library(arm)
> Loading required package: MASS
> Loading required package: Matrix
> [snip]
> Attaching package: 'arm'
>
> The following object(s) are masked from 'package:coda':
>
> traceplot
>
> > data("Mroz", package = "car")
> > fm <- glm(lfp ~ ., data = Mroz, family = binomial)
> > coefplot(fm)
>
> with version below.
>
> cheeres, Dave
>
> >
> > detach("package:arm")
>
> > coefplot <- function(object, df = NULL, level = 0.95, parm =
> NULL,
> + labels = TRUE, xlab = "Coefficient confidence intervals",
> ylab =
> "",
> + xlim = NULL, ylim = NULL,
> + las = 1, lwd = 1, lty = c(1, 2), pch = 19, col = 1,
> + length = 0, angle = 30, code = 3, ...)
> + {
> + cf <- coef(object)
> + se <- sqrt(diag(vcov(object)))
> + if(is.null(parm)) parm <- seq_along(cf)
> + if(is.numeric(parm) | is.logical(parm)) parm <-
> names(cf)[parm]
> + if(is.character(parm)) parm <- which(names(cf) %in% parm)
> + cf <- cf[parm]
> + se <- se[parm]
> + k <- length(cf)
> +
> + if(is.null(df)) {
> + df <- if(identical(class(object), "lm"))
> df.residual(object)
> else 0
> + }
> +
> + critval <- if(df > 0 & is.finite(df)) {
> + qt((1 - level)/2, df = df)
> + } else {
> + qnorm((1 - level)/2)
> + }
> + ci1 <- cf + critval * se
> + ci2 <- cf - critval * se
> +
> + lwd <- rep(lwd, length.out = 2)
> + lty <- rep(lty, length.out = 2)
> + pch <- rep(pch, length.out = k)
> + col <- rep(col, length.out = k)
> +
> + if(is.null(xlim)) xlim <- range(c(0, min(ci1), max(ci2)))
> + if(is.null(ylim)) ylim <- c(1 - 0.05 * k, 1.05 * k)
> +
> + if(isTRUE(labels)) labels <- names(cf)
> + if(identical(labels, FALSE)) labels <- ""
> + labels <- rep(labels, length.out = k)
> +
> + plot(0, 0, xlim = xlim, ylim = ylim, xlab = xlab, ylab =
> ylab,
> + axes = FALSE, type = "n", las = las, ...)
> + arrows(ci1, 1:k, ci2, 1:k, lty = lty[1], lwd = lwd[1], col
> = col,
> + length = length, angle = angle, code = code)
> + points(cf, 1:k, pch = pch, col = col)
> + abline(v = 0, lty = lty[2], lwd = lwd[2])
> + axis(1)
> + axis(2, at = 1:k, labels = labels, las = las)
> + box()
> + }
> >
> >
> > coefplot(fm, parm = -1)
>
>
>
>
> Achim Zeileis wrote:
>
> I've thought about adding a plot() method for the coeftest()
> function
> in
> the "lmtest" package. Essentially, it relies on a coef() and a
> vcov()
> method being available - and that a central limit theorem holds.
> For
> releasing it as a general function in the package the code is
> still
> too
> raw, but maybe it's useful for someone on the list. Hence, I've
> included
> it below.
>
> An example would be to visualize all coefficients except the
> intercept
> for
> the Mroz data:
>
> data("Mroz", package = "car")
> fm <- glm(lfp ~ ., data = Mroz, family = binomial)
> coefplot(fm, parm = -1)
>
> hth,
> Z
>
> coefplot <- function(object, df = NULL, level = 0.95, parm =
> NULL,
> labels = TRUE, xlab = "Coefficient confidence intervals", ylab
> = "",
> xlim = NULL, ylim = NULL,
> las = 1, lwd = 1, lty = c(1, 2), pch = 19, col = 1,
> length = 0, angle = 30, code = 3, ...)
> {
> cf <- coef(object)
> se <- sqrt(diag(vcov(object)))
> if(is.null(parm)) parm <- seq_along(cf)
> if(is.numeric(parm) | is.logical(parm)) parm <-
> names(cf)[parm]
> if(is.character(parm)) parm <- which(names(cf) %in% parm)
> cf <- cf[parm]
> se <- se[parm]
> k <- length(cf)
>
> if(is.null(df)) {
> df <- if(identical(class(object), "lm")) df.residual(object)
> else
> 0
> }
>
> critval <- if(df > 0 & is.finite(df)) {
> qt((1 - level)/2, df = df)
> } else {
> qnorm((1 - level)/2)
> }
> ci1 <- cf + critval * se
> ci2 <- cf - critval * se
>
> lwd <- rep(lwd, length.out = 2)
> lty <- rep(lty, length.out = 2)
> pch <- rep(pch, length.out = k)
> col <- rep(col, length.out = k)
>
> if(is.null(xlim)) xlim <- range(c(0, min(ci1), max(ci2)))
> if(is.null(ylim)) ylim <- c(1 - 0.05 * k, 1.05 * k)
>
> if(isTRUE(labels)) labels <- names(cf)
> if(identical(labels, FALSE)) labels <- ""
> labels <- rep(labels, length.out = k)
>
> plot(0, 0, xlim = xlim, ylim = ylim, xlab = xlab, ylab = ylab,
> axes = FALSE, type = "n", las = las, ...)
> arrows(ci1, 1:k, ci2, 1:k, lty = lty[1], lwd = lwd[1], col =
> col,
> length = length, angle = angle, code = code)
> points(cf, 1:k, pch = pch, col = col)
> abline(v = 0, lty = lty[2], lwd = lwd[2])
> axis(1)
> axis(2, at = 1:k, labels = labels, las = las)
> box()
> }
>
>
> On Fri, 2 Jul 2010, Tal Galili wrote:
>
> > Specifically this link:
> >
> http://tables2graphs.com/doku.php?id=04_regression_coefficients
> >
> > Great reference Bernd, thank you.
> >
> > Tal
> >
> >
> > ----------------Contact
> >
> Details:-------------------------------------------------------
> > Contact me: Tal.Galili at gmail.com | 972-52-7275845
> > Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il
> (Hebrew) |
> > www.r-statistics.com (English)
> >--------------------------------------------------------------------------
> -
> -------------------
> >
> >
> >
> >
> > On Fri, Jul 2, 2010 at 10:31 AM, Bernd Weiss <bernd.weiss at
> uni-koeln.de>wrote:
> >
> >> Am 02.07.2010 08:10, schrieb Wincent:
> >>> Dear all,
> >>>
> >>> I try to show a subset of coefficients in my presentation.
> It
> seems
> >>> that a "standard" table is not a good way to go. I found
> figure 9
> >>> (page 9) in this file (
> >>>
> >>http://www.destatis.de/jetspeed/portal/cms/Sites/destatis/Internet/DE/Con
> te
> nt/Wissenschaftsforum/Kolloquien/VisualisierungModellierung__Beitrag,proper
>
> ty=file.pdf
> >>>
> >>>
> >> ) looks pretty good. I wonder if there is any function for
> such
> plot?
> >>> Or any suggestion on how to present statistical models in a
> >>> presentation?
> >>
> >> Hi Wincent,
> >>
> >> I guess you are looking for "Using Graphs Instead of Tables
> in
> Political
> >> Science" by Kastellec/Leoni
> <http://tables2graphs.com/doku.php>.
> >>
> >> HTH,
> >>
> >> Bernd
> >>
> >> ______________________________________________
> >> R-help at r-project.org mailing list
> >> 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 at r-project.org mailing list
> > 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.
> >
>
> --
> Dave Atkins, PhD
> Research Associate Professor
> Department of Psychiatry and Behavioral Science
> University of Washington
> datkins at u.washington.edu
>
> Center for the Study of Health and Risk Behaviors (CSHRB)
>
>
> 1100 NE 45th Street, Suite 300
> Seattle, WA 98105
> 206-616-3879
> http://depts.washington.edu/cshrb/
> (Mon-Wed)
>
> Center for Healthcare Improvement, for Addictions, Mental
> Illness,
> Medically Vulnerable Populations (CHAMMP)
> 325 9th Avenue, 2HH-15
> Box 359911
> Seattle, WA 98104
> http://www.chammp.org
> (Thurs)
>
>
> ______________________________________________
> R-help at r-project.org mailing list
> 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.
>
>
>
>
>
>
More information about the R-help
mailing list