New answers tagged

1

In the two-stage approach, in the first step, you summarize the repeated measurements per peptide_Id into a single number. This inevitably leads to some information loss. To give another example, say that you measure the blood pressure of a patient ten times. Even though these measurements are correlated, they contain more information than their average. The ...


0

I'm assuming that s(time) or something like it in the model? If so, you can get into situations where both the smooth of time and the CAR(1) process, which are mathematically similar, cannot both be uniquely identified; one of the two ends up winning out, but which one can depend on a number of things. As you've found out it can depend on the scale of the ...


0

There are some pieces of an answer in the comments to the question, regarding the failure of hypothesis tests for the saturated model, and whether to handle some independent variables as factors. The following sketches out a potential solution in R for a post-hoc approach for a factorial anova. As an illustration, it uses different data for the BVTV ...


0

As I know, you need to have at least one replication in a $2^2 3$ mixed factorial design like yours. That's because if you don't have replications you not have degrees of freedom for the mean square of residuals. Maybe, you can collapse your design in two of the tree factor. Try with the moist factor.


1

1(a) Anova() can be easier to understand in terms of evaluating the significance of a predictor in your model, even though there is nothing wrong with the output from summary(). The usual R summary() function reports something that can appear quite different from Anova(). A summary() function typically reports whether the estimated value for each ...


0

One way to think about your question (or in general the philosophy regarding linear and non-linear models) is that all REAL models are non-linear. The only scenario where you will find a linear model is if you happen to simulate it yourself. When analyzing the data, especially after Econometrics 101, we start with the assumption that the model can be ...


0

You are fitting one and the same model using different starting parameters. The global optimum of a fit (assuming it exists) is independent of starting parameters. If you get different parameters, either these parameters combined result in the exact same fit or some of these fits did not converge to the global optimum. You could see this much better, if you ...


1

The lcmm() function fit a latent class linear mixed-effects model. This postulates that there are some underlying sub-populations in your data that you wish to recover. The model is estimated using maximum likelihood, and therefore it will provide you with correct inferences provides that any missing data in your outcome variable are missing at random and ...


0

The basic option, as in the effects package is to plot fitted values against one of the interacting regressors while fixing the other regressor at a representative value (say, mean, median or any "interesting" value). For effects package syntax and examples, see https://cran.r-project.org/web/packages/sjPlot/vignettes/plot_interactions.html Alternatively, ...


1

While the intuition behind it is not bad, I would not suggest using a Poisson glm for this problem, as it would have a different domain than your problem. You could find you glm predicting a value of 6 (which for you doesn't mean anything), and this can also wrongly affect the computation of the loss. I would instead suggest using a ordinal logistic (or ...


0

According to the documentation, you should get a the final parameter $\lambda$ from the train() function by accessing the bestTune value. In your case, that would be ridgeFit$bestTune. A fit object from the best parameter can be accessed with ridgeFit$finalModel. As an alternative, you should be able to get away with just calling summary(ridgeFit) and ...


1

Basically, the Wald statistic is not good and you shouldn't trust it for mixed models. It uses a much cruder approximation to the actual likelihood than you get with the profile and boot.ci methods. If R (and SAS and JMP and...) would have been written today, they would not have bothered implementing Wald stats. That's why the summary.merMod method ...


-1

All the details for the CP calculation with numerical examples are show here https://sites.google.com/stern.nyu.edu/rdeo/home


0

There are a number of feature selection techniques in random Forests. As Dij pointed out, RFE is a typical strategy used in random forests. Try the ones out below and see if that helps as well. These are just a few of the techniques available as R packages - some are easier to install than others, but give them a shot. AUCRF: Urrea, V., & Calle, M. (...


1

Normalizing the dependent variable as you have does not make sense. It would be one thing to take a regular z-score of it (based on the sample mean and standard deviation). Sometimes people do that to their predictors and outcome to get their coefficients into an effect size metric - 1 standard deviation increase in predictor is associated with XX standard ...


0

Your RMS method solution agrees with the solution that Frank Harrell approved in this answer. I've found survest to be extremely slow, particularly when running predictions for >10 bootstrap re-samples of test folds. But aside from code performance, the accuracy of the calculations should be more robust. HTH.


0

I would suggest to only transform your data when absolutely necessary. It's not best practice in the modern era when you're using a computer that can easily handle generalised linear models. Transformations inherently produce less useful data. Best thing to do is use the glm() function as you are, but to check the residuals from the different families of ...


0

I found a solution. Try s <- standardizedSolution(modelo.v1.fit) s[s$op == "~~",] to see the p values...


2

The package mcp was made just for scenarios like this. See below how I structured your data as df later. Fit a change point model First, let's define a slope followed by a joined plateau. We add varying (random) change point locations (the left-hand side of the equation): model = list( shoot ~ 1 + P, # intercept and slope 1 + (1|cultivar) ~ 0 # ...


2

Use mcp if you (1) want to quantify uncertainty about the location of the change point, and (2) want to specify a more informed model structure, e.g., that the first segment is a plateau. I arranged the data so that it is a regular data frame. Then fit an AR(1) model with a plateau + joined slope: model = list( y ~ 1 + ar(1), ~ 0 + x ) library(mcp) fit ...


0

For a Bayesian approach, you can use mcp to fit a Poisson or Binomial model (because you have counts from fixed-interval periods) with autoregression applied to the residuals (in the log space). Then compare a two-segment model to a one-segment model using cross-validation. Before we start, note that for this dataset, this model does not fit well and cross-...


1

The uncertainty in the mean is not just a rescaling of the total because there is uncertainty in the denominator as well. You get the mean by dividing the estimated total by the estimated population size, so there's a fairly complicated formula involving variances and covariances for the SE of the mean. A simple way to reproduce the calculation is to use ...


0

Appropriate for what purpose? That is the question. Since the residual deviance is significantly large, you can conclude that either your data is overdispersed or the model does not explain all the signal in the data. On the other hand, the model does explains a deviance of $165.8 - 29.92 = 135.9$ on 5 degrees of freedom, which is most ($135.9 / 165.8 = 82$...


0

I was looking for a solution that works for PCA performed using ade4. Please find the function below: library(ade4) irisX <- iris[,1:4] # Iris data ncomp <- 2 # With ade4 dudi_iris <- dudi.pca(irisX, scannf = FALSE, nf = ncomp) rotate_dudi.pca <- function(pca, ncomp = 2) { rawLoadings <- as.matrix(pca$c1[,1:ncomp]) %*% diag(sqrt(...


4

Great question! (Almost a professional statistician here!) Let's see if I can try to answer your question. Whuber gave a comment on what the typical strategy that someone might follow. I'll try to complement that by providing a couple papers that talk specifically about properties of tests when assumptions are violated. I'm doing this since, if I had a ...


0

You want to use the package betareg -> https://cran.r-project.org/web/packages/betareg/vignettes/betareg.pdf To install the package you can do: install.packages("betareg") library(betareg) And then like in a linear regression you could run: betareg(Perc_Reacting~Pulses+IndMutant, data) I will say tho that I never used this package or run a beta ...


0

One important point is what you mean by "a1:c4 at day 1 are different from a1:c4 at day 2 and so on". If this means that the a1 from day 1 is not more linked to a1 from day 2 that to a2 from day 2, then you don't have 10 simulations, but 10x30= 300 simulations (90 from category A, 90 from category B and 120 from category C). You should name them differently (...


0

A primary issue with your 468 daily values is the presence of a number of anomalies. I took your data into AUTOBOX ( which I have helped to develop) and it simultaneously idenifies both the arima component ( 1,0,0)(1,0,0)7 and a number of evidented anomalies. presents the Actual/Fit and Forecast. The model is here and here with statistics here The ...


3

It doesn't make sense to fit random slopes for Days without including it also as a fixed effect, unless you know the the overall effect of Days is zero. Try including it also a a fixed effect. If it still doesn't converge, try removing Days as a random slope, and just keeping it fixed. If the first measurement of Days is much above zero, then you might ...


2

By definition, 'white noise' is serially uncorrelated. It is like hitting piano keys at random; note n+1 cannot at all be predicted by note n. This is distinguished from brown noise and pink noise, both of which have some serial correlation, or dependence, albeit to different extents. To answer your question, 'strict white noise' assume essentially 0 serial ...


0

Your question isn't very clear about what exactly you wanted to do. That said, f_classif does ANOVA for feature selection for you instead of chi-squared. R package FSelector provides chi squared feature selection with function chi.squared. Here is the example from the document: https://rdrr.io/cran/FSelector/man/chi.squared.html


1

If you were to perform this as a linear regression, in R, you would receive the intercept term (the point at which the line touches the y-axis, which is equivalent to the value of y when the x-variable is 0). If you are using the ggplot2 package to graph the least squares line through the data points, ggplot can use a stat_smooth(method = "lm"), to help draw ...


3

Short answer: this is a code error. If I use the log scale, should I also use a log-proposal? The probability of acceptance is always $$\min\left\{1,\dfrac{\pi(x^\text{new})}{\pi(x^\text{old})}\times\dfrac{q(x^\text{old}|x^\text{new})}{q(x^\text{new}|x^\text{old})} \right\}$$ which can also be written as $$\min\left\{1,\exp[\ln\pi(x^\text{new})-\ln\...


3

Here is how to do this for one cultivar: plot(shoot ~ P, data = subset(DF, cultivar == "Dinninup")) fit1 <- nls(shoot ~ ifelse(P < bp, m * P + c, m * bp + c), data = subset(DF, cultivar == "Dinninup"), start = list(c = 1, m = 0.05, bp = 25), na.action = na.omit) summary(fit1) #Formula: shoot ~ ifelse(P < bp, m * P + c, m * ...


2

I found this question after a couple of months, so I assume a solution was already found. Nevertheless, here an example, just for the record: library(deSolve) ## matrix multiplication makes this model compact multi <- function(t, n, parms) { with(parms, { dn <- r * n + n * (A %*% n) return(list(c(dn))) }) } times <- seq(from = 0, to = ...


0

I'll answer your questions in reverse. Your model posits that the log odds of "jumping the hurdle" look like $$ \operatorname{logit}(p) = \hat{\beta}_0 + \hat{\beta}_1\operatorname{Start}$$ With $\hat{\beta}_0 = 0.51$ and $\hat{\beta}_1 \approx 0$. To obtain the binomial probability, invert the logit using $$ p = \dfrac{1}{1+\exp(-(\hat{\beta}_0 + \hat{\...


0

It depends on how you build your prediction model. For example, suppose you have a model that says: for region A the score is 1 and for all other regions the score is zero. Then your model will correctly predict the score for the first two people and incorrectly for the others. Or you can have a model that predicts for region A a score of 0.95, for region B ...


2

As a comment - note that the Cauchy distribution does not have a mean, so I interpret your reference to the mean as being to the location parameter of the Cauchy, which, due to the symmetry of the Cauchy, is also the median of the distribution. With respect to bias: the sample median is an unbiased estimator of the population median, so you should not ...


3

Sounds like regularization models could be of use to you (Lasso especially, since your problem seems to be centered around picking variables, as opposed to Ridge and Elastic net). However, these can be difficult to set up. They are included in Stata 15 if you have access to that. In R, it's often a bit more handywork but hopefully someone can chime in with ...


1

I am assuming (perhaps wrongly) that intensity is recorded in some standard way (e.g. like weight being recorded in kilograms). In this case, the measures are on the same scale now. If you standardize, they will be on different scales, but each sample will have the same mean and sd. Which do you want? That depends on what you are trying to do. Why do the ...


1

There are two types of Intervention studies . The first one is called Intervention Analysis (de jure )..the second is called Intervention Detection (de facto). Simply search here for R and one or the other. The ultimate approach is to use a SARMAX model https://autobox.com/pdfs/SARMAX.pdf to form a useful equation leading directly to tests of statistical ...


17

As mentioned by @DemetriPananos, theoretical justification would be the best approach, especially if your goal is inference. That is, with expert knowledge of the actual data generation process, you can look at the causal paths between the variables and from there you can select the variables which are important, which are confounders and which are mediators....


4

Theoretical justification beyond all else. Aside from that, LASSO or similar penalized methods would be my next suggestion.


0

Just wanted to add that you can do this with many other packages. If you want to get an estimate of uncertainty around the change point (something nls cannot do), try the mcp package. # Simulate the data df = data.frame(x = 1:100) df$y = c(rnorm(20, 50, 5), rnorm(80, 50 + 1.5*(df$x[21:100] - 20), 5)) # Fit the model model = list( y ~ 1, # Intercept ~ ...


0

Depending on what you are trying to estimate, sometimes it is useful with imbalanced data to re-sample the data. There are various algorithms out to resample your data. You can undersample your zeros to prevent them from dominating the model. You can oversample (aka duplicate) your non-zero values to make them have a higher impact on your model. You can use ...


5

I assume x is a predictor, not the dependent variable. Fit x as quadratic with an exception that allows for a discontinuity at zero, i.e., add an indicator variable for x > 0.


0

I think the answer is triangulation. Triangulate all trees and then break the edges where one vertex is dead tree and another vertex is surviving tree.


11

As a point for further elaboration, here is the explanation in the thread you reference: If it's mathematically possible, [this method] will find an $X_{Y_1,Y_2,\ldots,Y_k;\rho_1,\rho_2,\ldots,\rho_k}$ having specified correlations [between $X$ and] an entire set of $Y_i$. Just use ordinary least squares to take out the effects of all the $Y_i$ from $X$ ...


0

R's forecast package now has a function mstl() to handle multiple seasonal time series decomposition. This page has got more details how to use it: https://pkg.robjhyndman.com/forecast/reference/mstl.html


1

As a supplement to whuber's answer, I've written a script in Python which goes through each step of the sampling scheme. Note that this is meant for illustrative purposes and is not necessarily performant. Example output: n=10, s=20, k=4 Starting grid . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . X X X X X X X X X X Filled in grid X X . . X ...


Top 50 recent answers are included