Re: [R] How does predict() calculate prediction intervals?

2013-01-30 Thread Cade, Brian
Kurt: You missed including the term 1/n in your var.y.hat calculation. Do that and pred and hand are the same. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: brian_c...@usgs.gov tel: 970 226-9326

Re: [R] When the interaction term should be interpreted in AIC table?

2013-05-24 Thread Cade, Brian
Galina: The AIC, delta AIC, and AIC weights all reference an entire model and provide no information on how you should interpret the individual parameters within a model. If you believe based on your AIC weights that the model with ZONE + LABRADOR TEA + YEAR + ZONE x LABRADOR TEA is a

[R] indexing names for looping across computations done on pairs of matrices

2014-04-15 Thread Cade, Brian
So I know I must be missing something simple and obvious for the following data manipulation where I have (in this example) 11 pairs of matrices (gs4.0 to gs4.100 and ps1.0 to ps1.100) from some population simulations (all with same dimensions) where I want to get some summary statistics on the

Re: [R] indexing names for looping across computations done on pairs of matrices

2014-04-16 Thread Cade, Brian
, April 15, 2014 4:30 PM, Cade, Brian ca...@usgs.gov wrote: So I know I must be missing something simple and obvious for the following data manipulation where I have (in this example) 11 pairs of matrices (gs4.0 to gs4.100 and ps1.0 to ps1.100) from some population simulations (all with same

Re: [R] Nonparametric MANOVA

2014-05-13 Thread Cade, Brian
Rich: Your use of the term MANOVA suggests a multivariate response (Y). If what you really have is multiple factors (predictors), then this is a different modeling construct (multiple regression) and it would seem nonpartest() is not appropriate. I've been analyzing water quality constituents

Re: [R] Interpreting results of rq() when tau = -1

2014-05-14 Thread Cade, Brian
Rich: When you use tau=-1 in rq() then under the $sol you are getting estimates for all possible quantiles for the given regression model, which includes the objective function value minimized, the quantile value, estimates for all the parameters, and I think the output still includes the

Re: [R] Problem with zero-inflated negative binomial model in sediment river dynamics

2013-08-13 Thread Cade, Brian
Lauria: For historical reasons the logistic regression (binomial with logit link) model portion of a zero-inflated count model is usually structured to predict the probability of the 0 counts rather than the nonzero (=1) counts so the coefficients will be the negative of what you expect based on

Re: [R] Problem with zero-inflated negative binomial model in sediment river dynamics

2013-08-14 Thread Cade, Brian
Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326 On Wed, Aug 14, 2013 at 4:07 AM, Achim Zeileis achim.zeil...@uibk.ac.atwrote: On Tue, 13 Aug 2013, Cade, Brian wrote: Lauria: For historical reasons the logistic regression (binomial

Re: [R] Nonparametric k-way ANOVA

2013-10-24 Thread Cade, Brian
Depends on what the original poster really wants with a nonparametric k-way ANOVA. If what they are really looking for is an approach that has fewer distributional assumptions, then it is possible to perform permutation tests for linear models parameterized for k-way anova that would eliminated

Re: [R] Beta-coefficients for ZINB model

2012-12-14 Thread Cade, Brian
I always liked Johan Bring's (1994. How to standardize regression coefficients. Am. Stat. 48(3):209-213.) arguments for standardizing the beta estimates in a linear model by their partial standard deviations which equates to comparing ratios of t-statistics between variables to determine their

Re: [R] NADA and cenmle

2013-03-14 Thread Cade, Brian
Shane: Just to add some practical advice on top of Rich's interpretation of the censoring process (which is correct), my recent experiences with analyzing below-detection limit chemical concentrations in water using left-censoring estimators indicates that historically people have not always

Re: [R] How many samples ACTUALLY used in regression?

2013-03-18 Thread Cade, Brian
Perhaps a crude but reliable way is to check the number of residuals, e.g., length(my.model$resid). Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326

Re: [R] Writing contrast statements to test difference of slope in linear regressions

2013-04-23 Thread Cade, Brian
Tyler: To compare estimates of different slopes for Time by different factor levels of Temperature you need to estimate a model that allows for separate slopes. What you've estimated is a model for just different intercepts associated with Temperature factor but a common slope for Time. You

Re: [R] MuMIn Random Effects Variance

2013-12-13 Thread Cade, Brian
You might want to rethink about getting model averaged coefficients. That is a bunch of nonsense if you have any multicollinearity among the predictors. Model averaged predictions might be useful. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave.,

Re: [R] Simple permutation question

2014-06-25 Thread Cade, Brian
It is called sample(,replace=F), where the default argument is sampling without replacement. Try x - c(A,B,C,D,E) sample(x) Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov

Re: [R] Percentage cover data with many zeros

2015-01-26 Thread Cade, Brian
Ben: You have a statistical problem with a bounded response variable (0 to 100%, or 0.0 to 1.0) and thus, might make use of a logistic quantile regression model (see Bottai et al. 2010. Logistic quantile regression for bounded outcomes. Statistics in Medicine 29: 309-317). This requires a

Re: [R] Drawing the regression line and the 95% confidence intervals

2015-05-06 Thread Cade, Brian
The prediction intervals are likely to be much wider than the confidence intervals so you will need to be sure you scale the yaxis limits large enough to see them. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818

Re: [R] alternatives to KS test applicable to K-samples

2015-06-01 Thread Cade, Brian
email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326 On Sat, May 30, 2015 at 3:11 PM, David Winsemius dwinsem...@comcast.net wrote: On May 29, 2015, at 10:02 AM, Cade, Brian wrote: Wensui: There are the multi-response permutation procedures (MRPP) that readily test the omnibus

Re: [R] alternatives to KS test applicable to K-samples

2015-05-29 Thread Cade, Brian
Wensui: There are the multi-response permutation procedures (MRPP) that readily test the omnibus hypothesis of no distributional differences among multiple samples for univariate or multivariate responses. There also are empirical coverage tests that test a similar hypothesis among multiple

Re: [R] monte carlo simulations in permanova in vegan package

2015-10-27 Thread Cade, Brian
Sean: There are only 20 possible combinations, 6!/(3! x 3!), so you just need to enumerate them completely (no Monte Carlo approximation required). I don't know if permanova() can do this but you can do it with the mrpp() functions and argument (,exact=TRUE) in Blossom package for R. Brian

Re: [R] algorithmic method quantile regression

2015-10-15 Thread Cade, Brian
>From ?rq.fit.pfn I see: Details: Preprocessing algorithm to reduce the effective sample size for QR problems with (plausibly) iid samples. The preprocessing relies on subsampling of the original data, so situations in which the observations are not plausibly iid, are likely

Re: [R] Comparing two populations based on the percentile values calculated from two independent samples

2015-09-24 Thread Cade, Brian
You can use quantile regression (in quantreg package) to compare percentiles between two groups in a linear model formulation. If you are really interested in "equivalence" testing and not just using the term informally, you might check out Cade (2011. Estimating equivalence with quantile

Re: [R] Right censored data, abundant in zeros for regression analysis.

2015-12-28 Thread Cade, Brian
Tom: One possibility might be to use the censored quantile regression implementation (crq) in the quantreg package (accommodates left or right censoring) across a range of quantiles (e.g., 0.05 to 0.95) but where interest is likely to be focused on estimates for quantiles greater than the

Re: [R] Power Calculation:2-sided exact equivalence test for Binomial Proportions

2016-06-07 Thread Cade, Brian
MJ: I think the EnvStats package has various power functions for binomial applications (also confidence interval half-widths). Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov

Re: [R] lm() with spearman corr option ?

2016-04-29 Thread Cade, Brian
I think you would just need to replace the lm() function call with cor(x,y,method="spearman". It would probably be more informative to actually plot by the magnitude of the correlation coefficient (all |r| >= 0.20 or something similar) rather than just by those with P <=0.05. Brian Brian S.

Re: [R] Heteroscedasticity in a percent-cover dataset

2016-04-15 Thread Cade, Brian
Quit trying to eliminate heteroscedasticity in your data - there is information there in the pattern of changing variances. I would suggest instead that you go directly after modeling the change in entire distributional form by using quantile regression (package quantreg). So, for example,

Re: [R] Linear model vs Mixed model

2016-07-13 Thread Cade, Brian
ict(fit_lm2)) >123456 > 295.0310 305.4983 315.9656 326.4329 336.9002 347.3675 > > And these are quite different from the mixedmodel: > > > fit_lmer = lmer(Reaction ~ Days + (1| Subject), sleepstudy) > > head(predict(fit_lm

Re: [R] glmmLasso with interactions errors

2016-07-14 Thread Cade, Brian
It has never been obvious to me that the lasso approach can handle interactions among predictor variables well at all. I'ld be curious to see what others think and what you learn. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort

Re: [R] Linear model vs Mixed model

2016-07-12 Thread Cade, Brian
Your lm() estimates are using the default contrasts of contr.treatment, providing an intercept corresponding to your subject 308 and the other subject* estimates are differences from subject 308 intercept. You could have specified this with contrasts as contr.sum and the estimates would be more

Re: [R] non-parametric manova with post-hoc test

2017-01-18 Thread Cade, Brian
You could try a multi-response permutation procedure (MRPP) for multivariate hypothesis testing (null is groups come from a common distribution) without resorting to ranks. There are no automated multiple comparison procedures, but one could either look at pairwise contrasts of group (if that is

Re: [R] R-square prob is not calculated by randomization in lmPerm::lmp

2016-09-06 Thread Cade, Brian
For a linear model without an intercept term as in this example, neither the usual permutation scheme for testing Ho: B1 = 0 nor usual definition of R-squared apply. So you need to check what the developer of this code chose to do. If I'm recalling correctly, in a linear model with an intercept

[R] box type in Hmisc xYplot

2016-09-15 Thread Cade, Brian
Does anyone know how to change the box type in Hmisc package function xYplot. I want only the left and bottom axes drawn, similar to what I would accomplish with bty="l" argument in plot() function. bty= argument did not do anything for me in xYplot(). Brian Brian S. Cade, PhD U. S.

Re: [R] box type in Hmisc xYplot

2016-09-15 Thread Cade, Brian
coming along > and sticking things into it." > -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip ) > > > On Thu, Sep 15, 2016 at 1:32 PM, David Winsemius <dwinsem...@comcast.net> > wrote: > > > >> On Sep 15, 2016, at 12:59 PM, Cade, Br

[R] transposing x and y axes in xYplot

2016-08-23 Thread Cade, Brian
Is there a simple way to transpose the x and y axes with the xYplot() function in the Hmisc package, where y is a vector of point estimate and lower and upper confidence interval endpoints? What I'm looking for is something akin to coord_flip() used with ggplot(). Brian Brian S. Cade, PhD U.

Re: [R] transposing x and y axes in xYplot

2016-08-23 Thread Cade, Brian
Cheers, > Bert > > > Bert Gunter > > "The trouble with having an open mind is that people keep coming along > and sticking things into it." > -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip ) > > > On Tue, Aug 23, 2016 at 2:42 PM, Cade,

Re: [R] Confidence interval on quantile regression (quantreg)

2016-10-19 Thread Cade, Brian
Depending on the procedure used for estimating the CI, especially if the default rankscore inversion method, then it is possible that legitimate end points of the intervals for some quantiles with a given alpha (e.g., 0.05 for 95% CI) cannot be refined beyond plus or minus infinity. Of course,

Re: [R] GAM with the negative binomial distribution: why do predictions no match with original values?

2016-11-22 Thread Cade, Brian
Well part of the issue is that the negative binomial estimates are for means and they can differ a fair bit from the raw counts, but I'm also guessing that part of the issue is that the offset may not be accounted for with the predict.gam() function. Brian Brian S. Cade, PhD U. S. Geological

[R] warnings about factor levels dropped from predict.glm

2017-12-05 Thread Cade, Brian
I am helping a student with some logistic regression analyses and we are getting some strange inconsistencies regarding a warning about factor levels being dropped when running predict.glm(, newdata = ournewdata) on the logistic regression model object. We have checked multiple times that the

[R] glm$effects

2018-01-12 Thread Cade, Brian
I know I must be missing something obvious, but checking help and googling a bit did not turn up a useable answer. When I've estimated a glm() model object (my example is with just identity link with gaussian family so I could have used lm() instead), one of the terms returned in the model object

Re: [R] question regarding the AICcmodavg package

2018-02-20 Thread Cade, Brian
among predictors and some don't. Stick to model averaging the predicted responses and you can do something that is sensible and that can be applied to any combination of predictor variables including those with interactions. Brian Cade Brian S. Cade, PhD U. S. Geological Survey Fort Collins

[R] overlaying graphs in xYplot (Hmisc)

2020-04-22 Thread Cade, Brian S via R-help
Hi All. I am trying to construct a graph using the xYplot() function in Hmisc package (thank you Frank Harrell) taking advantage of the Cbind() argument for plotting the median, 10th, and 90th quantiles and also the cbind() argument for individual data values. I know how to do both of these

Re: [R] [EXTERNAL] Re: overlaying graphs in xYplot (Hmisc)

2020-04-22 Thread Cade, Brian S via R-help
Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov<mailto:brian_c...@usgs.gov> tel: 970 226-9326 From: David Winsemius Sent: Wednesday, April 22, 2020 11:10 AM To: Cade, Brian S ; r-help Subject: [EXTERNAL] Re: [R] overlaying graphs in

[R] complicated time series filtering issue

2022-04-04 Thread Cade, Brian S via R-help
Hello: I have an issue with filtering in a time series of animal growth data that seems conceptually simple but I have not come up with effective code to implement this. I have temporal sequences of lengths by individuals and I want to retain only those data that are >10 days apart

Re: [R] [EXTERNAL] Is there a canonical way to pronounce CRAN?

2022-05-04 Thread Cade, Brian S via R-help
CRAN as in cranberry. Just my minor thought. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov tel: 970 226-9326 From: R-help