title: Problem in understanding hidden Markov model EM training and log-lik computation result

text: i'm using Kevin Murphy's HMM Toolbox. But i have a problem in understanding the results. My problem is to classify some sequences with Hidden Markov Model. I'm using a (test) dataset build with random sequences of symbol '1' '3' '6'. this one is my training dataset: data = [ 1 1 1 3 6; 1 1 3 3 6; 1 3 3 6 6] - 3 observations of 5 symbols. i use EM baum welch to train and log_prob to evaluate this new observation: data1 = [1 3 6] i'm expecting logprob near 0 like -0.1 instead i have thoose result:               loglik = -110.04     loglik = -110.63     loglik = -91.679      why have i those weird result?? this is my MATLAB code:               O = 10;     Q = 12;          % training data     T = 5;     nex = 3;     data = [ 1 1 1 3 6;              1 1 3 3 6;              1 3 3 6 6]          % initial guess of parameters     prior1 = normalise(rand(Q,1));     transmat1 = mk_stochastic(rand(Q,Q));     obsmat1 = mk_stochastic(rand(Q,O));          % improve guess of parameters using EM     [LL, prior2, transmat2, obsmat2] = dhmm_em(data, prior1, transmat1, obsmat1, 'max_iter', 50);     LL          % use model to compute log likelihood     data1 = [1 3 6]     loglik = dhmm_logprob(data1, prior2, transmat2, obsmat2)      EDIT: also giving the same data of testing i recive weird result! why that?


title: Problems with interpretation of a CLASS variable with many levels in a logistic regression model

text: Thank you in advance. I am trying to build a model to predict the possibility of ordinal accident levels in a nuclear station if an accident happens. The accidents include four levels: Site_Area_Emergency, Alert, Unusual Event and Non_Emergency (In order of from Strongest to Weakest, Non Emergency is the weakest). I have three variable as potential predictor variables: The Age of that power station and the Facility_or_Organization type of that power station and power output. The age variable and power output are numerical variable. The Facility_or_Organization is obviously a nominal variable (like Kansas, New York, etc.) Therefore, based on the stuff I described above, I chose to build a Ordinal Logistic Model. I use the backward selection method. It turns out that only age and Facility_or_Organization type matter. My report is attached as a pdf file. It is a pdf file and quite long, therefore, I uploaded it to a share website, here is the link. Sorry for disturbing. http://minus.com/mvsqZH/43f My question is that: the overall degree of freedom is 67 because we have 65 different Facility_or_Organization types. The probability > chi_square for facility overall is quite small like 0.0001 but for each Facility_or_Organization type it is quite high. I took screen shots for you to have a look. ![enter image description here](http://i.stack.imgur.com/TP6JU.jpg) **My Question:** May be you already noticed that in the second last table called "Type 3 analysis of Effects". The Pr>ChiSq for Facility_Or_Organization is <0.001. However, in the last table called "Analysis of Maximum Likelihood Estimate" the Pr>ChiSq for each Fasicility is quite high around 0.96. For example, for first Facility ARKANSAS NUCLEAR, the Pr>ChiSq is 0.9670. I just want to know why the Pr>ChiSq for whole Facility_Or_Organization is quite low but for each of them it is high?


title: Random intercept mixed model in R - is this fit... OK?

text: I am looking for opinions/interpretation on a model I'm trying to fit. A disclaimer: I'm an ecologist and not a statistician; and I appreciate everyone's time and opinions! I am generally following Zuur et al. 2009's (using R) recommendations in an attempt to model travel speed as a function of a few environmental covariates, the main one of interest being a 1/0 "treatment" variable. The data are by nature are hierarchical and I have therefore used a mixed model with a random intercept. The mixed model fits significantly better than the base `gls` model. Additionally I included a `varIdent` with the mixed model, and that is significantly better than without. From this point, I have looked at assumptions and diagnostics - the residuals look OK, I think, but I always have a hard time telling if it's "too much" pattern. I think less than 5% of the data are outside of the 2/-2, so I'm good there. The variation doesn't seem to be driven by any points in particular, and the variance of residuals for each group of the "treatment" is close to 1 (using `tapply`; I interpret this as minor heterogeneity). A box plot of the residuals by the random effect shows equal spread around zero and minor variation = independence (see below)? I've looked at residuals vs. fitted (not too bad), residuals vs. each covariate (pattern with one of my covariates, elevation), histogram of the residuals (still quite skewed), and have checked for all potential influential observations. I just don't know if I'm "done" and can now consider this model an OK fit? I've included a couple of diagnostic images (response variable can't go below a certain value, and therefore the residual plot shows a sharp "line" in the negative values) - std residuals vs fitted; residuals by random intercept; histogram of residuals (not normalized); residuals by elevation covariate - could be something going on here, but not sure how to deal with it: ![enter image description here](http://i.stack.imgur.com/1CzXF.png) ![residuals by random intercept](http://i.stack.imgur.com/w3wJh.png) ![histogram of residuals \(not normalized\)](http://i.stack.imgur.com/G9XBz.png) ![resids by elevation covariate](http://i.stack.imgur.com/IlWVS.png)


title: How can I find out if a subset of Stack Exchange users increase/decrease their post rate based on badges earned?

text: I'm trying to mine that Stack Exchange data dump to find out whether there is a cluster of users that may be positively or negatively affected by the number of badges they've been awarded. The theory I'm working on is that for some people who were motivated by badges, that they eventually lose their appeal and they stop posting. Currently, I've filtered the raw data down to a set of:   * User ID   * Year and Month   * Posts made during that month   * Total badges that user has earned up to that month I'm not sure how to proceed in order to find these users. I'm not even sure how to formulate this question into something I can take action on. If I visualize it in my head, what I'm thinking is seeing a chart of an x-axis of badges earned, a y-axis of posts made per month, and then seeing the line dip after a certain number of badges. But this only works on a per user level. I want to find the cluster of users for whom that chart would hold. Any ideas? Thank you!


title: $\chi^2$ test - correcting for not fully independent sample

text: I have carried out a study in which I have gathered 3 sentences per participant. These sentences were then classified in 2 ways. I want to test if there are significant interactions between the classifications and find out the effect size. I have done a chi square test for independence and calculated Cramer's $V$ (on the counts of sentences for the cross of the 2 classifications). However, I am concerned that this is not correct, since the sentences are not fully independent (because each participant has entered 3 sentences). Is there a good way to correct for this (or are they other statistical tests that I should carry out instead of $\chi^2$ and/or Cramer's $V$)? A colleague of mine suggested removing the sentences that were entered by participants who had sentences in more than 1 category of a classification.


title: Time series with multiple subjects and multiple variables

text: I am a web developer and novice statistician. My data looks something like this               Subject  Week   x1  x2  x3  x4  x5  y1     A        1      .5  .6  .7  .8  .7  10     B        1      .3  .6  .2  .1  .3  8     C        1      .3  .1  .2  .3  .2  6       A        2      .1  .9  1.5 .8  .7  5     B        2      .3  .6  .3  .1  .3  2     D        2      .3  .1  .4  .3  .5  10        I am trying to predict y1 as a product of the x variables. However, I have reason to believe that there may be a lag in the effect of the multiple x variables on y1, i.e the x variables from week 1 for subject A influence y1 for subject A in week 2. Note that not all subjects will have data points for every week (in fact most won't). Subjects will tend to have data points for say week 1, 2, 3, 4 then drop off and not show up again until week 7,8,9. I am willing to restrict my analysis to data points where we have data for the previous N weeks given my hypothesis about lag. Like I said, I am a novice and am unsure of the best way to deal with a dataset of this form. I am hoping to carry out this analysis either in R, Python, or some combination of the two. I don't think that the current week's x variables will have no effect. I think they will have some effect, perhaps greater than previous weeks. I just believe that previous weeks will have some effect. I am expecting there to be two to three weeks of lag. To give a little context, the analysis that I am attempting here relates to judging the quality of online traffic. Every week I get a score rating the quality of a certain stream of users I send to a given website. I am trying to find secondary metrics, such as browser distribution, percent duplicate clickouts, etc. that will allow me to predict what that score will be ahead of time.


title: Estimating Frechet distribution parameter in R

text: I have a data set with 925 members that using Easyfit, I've determined is best fit using a Frechet distribution. Now, I'm trying to replicate this functionality in R to determine the optimum fit parameters using the VGAM package. Using this code:               vglm(V1 ~1, frechet2(location = 0), data, trace = FALSE)      I get the following error:               Error in eval(expr, envir, enclos) : initial values for 'location' are out of range      As I'm trying to estimate the 2-parameter distribution the location should be 0 (correct?) In the VGAM documentation it mentions that the package may fail for low values of shape parameter (near 2 or lower); however, in EasyFit the estimation shows the shape parameter at .2646. Any thoughts on how I should handle this? I've done a search on but haven't found any other packages that provide the same functionality. Any help would be welcome!


title: Useful heuristic for inferring multicollinearity from high standard errors

text: I was recently looking at some standard error calculations using the following:   * the formula for the covariance matrix in ordinary least squares (OLS)   * the "robust" formula for the heteroskedasticity-consistent covariance matrix   * the "modified-robust" formula (also called 'White' for its discoverer) for a bias-corrected and heteroskedasticity-consistent covariance matrix For a coefficient estimate that is about -0.05, the OLS standard error is about 0.0198. Meanwhile, the other two are in good agreement with each other and both are basically equal to 1.831. When I re-run the calculations using different libraries (SciPy.stats, scikits.statsmodels, and in R) I get excellent agreement, so I believe these estimates are correct. This feels like a really big difference to me, relative to the value of the coefficient, but I don't really know why it feels like a big difference and don't have a good quantitative way for interpreting the difference in standard errors. I know that if there is multicollinearity, it often means there are high standard errors. Under what conditions (if any) is it valid to go in the other direction and claim that if there is a high standard error, then there is probably multicollinearity? In my case, the regression is a very simple two-term linear model with a constant and a coefficient on an indicator variable. But there definitely could be correlation between the constant (a base wage level) and the indicator (veteran status). **Full disclosure:** this is related to a homework problem, but the homework problem just asked us to calculate the above standard errors. It didn't ask us to provide any interpretation, and I am just curious. If this had just been a data set from the wild and I had seen such differences in standard errors, it probably would make me feel nervous that I'm doing something wrong. I'm asking here just to see what a principled way to actually interpret the numbers might be.


title: Likes/interests datasets for social networks

text: I know there are publicly available datasets for social network graphs (I'm aware of an existing question for that). But does anyone know if there are any datasets available for user "likes" or interests (Facebook or any other social network)?


title: Criticism of Pearl's theory of causality

text: In the year 2000, Judea Pearl published _Causality_. What controversies surround this work? What are its major criticisms?


title: Parallel Stochastic Gradient Descent

text: I am trying to solve a linear regression problem with multiple variables in a Hadoop cluster. I was thinking of using SGD to estimate the parameters, but there is one issue I don't know how to solve. The SGD is calculated in parallel through all the nodes of the cluster, and my question is how to combine the output of those nodes. I mean, if I have for example a data set of 1000 rows and 5 nodes estimating the parameters using SGD with 200 rows of the data set each of them, how would I combine in a last step those 5 outputs?


title: uncorrelated random effects with lmer and glmer

text: I am modelling the activity of animals (subjects). Each subject is measured for activity repeatedly (which requires a random intercept associated with subject). Subjects have also been measured in two environments (a and b). I specifically want to see if the activity level of a given subject in environment a correlates with its activity level in environment b. In line with this really useful answer (How to specify uncorrelated random slopes in lmer() syntax?), I tried doing a log-likelyhood ratio test on these models:               full.model <-          lmer(activity ~ (1+environment|subject), family='binomial')     no.correlation.model <-          lmer(activity ~ (1|subject) + (environment-1|subject), family='binomial')      It appears that this code works nicely for gaussian data, but that it does not for binomial data. In this case, no.correlation.model has two random intercepts, a 'slope' associated with environment and a correlation. Is there a way to test whether the correlation estimated in full.model is significant? Here is a reproducible dataset.               activity <-rbinom(n=1000, size=1, prob=0.5)     subject <- rep(1:100,each=10)     environment <- rep(c('a', 'b'), 500)      Many thanks in advance!


title: What are my options for estimating probabilities of success/failure?

text: So, I want to estimate the probability that a student will succeed or fail in a particular major (suppose success means in this case to graduate with that major) based on data from that student's application. I have a lot of records of past students' application data and success. I am aware that I could try logistic regression, but I would like to know all my options, including methods which allow more easily for non-linearities. If I were just trying to classify incoming students, I think I would try an SVM model. But I don't know how easy it is to adapt that method to the task of estimating probabilities. ### EDIT: Much of the data is categorical. Often binary. Some features are real numbers / integers, as well.


title: What is the dimension (k) of these regression models?

text: I am attempting to use Akaike's Information Criterion to select the most appropriate model for some data. This means I need to find the likelihood of my data under various models, compute the AIC score for each model, and then select the model with the lowest AIC score. To find AIC, I would use: $$ AIC = 2k−2\log(L) $$ I want my question to be as clear as possible, so I'm going to go step-by step in calculating the likelihood of the data, simply by quoting from Barlow (1993). The data are a set of (x,y) pairs, exactly as on p. 93 of Barlow (1993), $\\{(x_i,y_i)\\}$. "The $x_i$ are known exactly, whereas the $y_i$ have been measured, each with some known resolution $\sigma_i$. Suppose that $y$ is given by some function $f$ of $x$, which also depends on a parameter (or, in general, several parameters) $a$, and that $a$ is what we want to estimate. This is the ideal $y$, however, and the actual measurements have been smeared by the resolution. "Invoking the central limit theorem to show that the distributions of the measured $y$ values about their ideals is Gaussian, the probability of a particular $y_i$, for a given $x_i$, is "$P(y_i;a) = \frac{1}{\sigma_i\sqrt{2\pi}}e^{-[y_i-f(x_i;a)]^{2}/2\sigma_i^{2}}$ "The logarithm of the likelihood for the complete data set is thus "$\ln{L}=-\frac{1}{2}\sum[\frac{y_i-f(x_i)}{\sigma_i}]^2-\sum\ln{\sigma_i\sqrt{2\pi}}$ " I have three models to select among: $$ y_A = m_1x $$ $$ y_B = m_1x^{2} $$ $$ y_C = m_1x^{m_2} $$ where $m_i$ is a parameter that I estimate using Least Squares minimization in KaleidaGraph. I'm clear on how to calculate the likelihood $L$ of the data under each of these three models, following the parameter estimation. Based on responses to an earlier version of this question, there may be some disagreement/confusion over how I account for measurement error. I have, I think, assumed that $s_i = \sigma_i$. If this is wrong, please explain. I am unsure of which values to use for $k$. For $y_A$ and $y_B$ I have estimated only one parameter, so does this mean that $k_A=1$ and $k_B=1$ ? Does $k_C=2$ ? Suppose the models included lower order terms as well, eg $$ y_A = m_1 + m_2x $$ $$ y_B = m_1 + m_2x + m_3x^{2} $$ $$ y_C = m_1 + m_2x^{m_3} $$ When selecting from this set, would I use $$ k_A = 2 $$ $$ k_B = 3 $$ $$ k_C = 3 $$ References: Barlow, R., 1993. Statistics: A Guide to the use of Statistical Methods in the Physical Sciences


title: Optimal number of components in a Gaussian mixture

text: So, getting an "idea" of the optimal number of clusters in k-means is well documented. I found an article on doing this in gaussian mixtures, but not sure I am convinced by it, don't understand it very well. Is there a ... gentler way of doing this? Greetings


title: Looking for a test for shape comparison

text: I have two different time serie both length = 100 and I need to know what is the best test (non-parametric, if possible), that return how much these two series are same or similar shapes. Here are two examples: first one the two series are very similar shape, not matching perfectly but overall trends and tops and bottoms match. ![enter image description here](http://i.stack.imgur.com/4GD6D.png) Second example, these two serie are different both for tops, bottoms and trends ![enter image description here](http://i.stack.imgur.com/BoPqH.png) Thanks all **Edit:** I am testing my series with Spearman's correlation coefficient, but I am not sure if this test can suits my need since Spearman test use a monotonic function, while my series are non-monotonic. Am I correct?


title: What optimization problem does least angle regression try to solve?

text: In Hastie et al's Elements of Statistical Learning, it says > Least angle regression (LAR) ... can be viewed as a kind of “democratic” > version of forward stepwise regression (Section 3.3.2). As we will see, LAR > is intimately connected with the lasso, and in fact provides an extremely > eﬃcient algorithm for computing the entire lasso path as in Figure 3.10. Forward stepwise regression tries to solve the following optimization problem for selection of best subset of features: $$ \min_x \|Ax-b\|_2 $$ s.t. $$ \|x\|_0 \leq M. $$ LASSO tries to solve the following optimization problem $$ \min_x \|Ax-b\|_2 $$ s.t. $$ \|x\|_1 \leq \epsilon. $$ I was wondering what optimization problem does least angle regression try to solve? Thanks!


title: Is a bounded normal distribution the same as a poisson distribution (for continuous variables)?

text: Lets say I have a normal distribution of estimates (e) for a random variable (v). At some point an event occurs where the value of the variable becomes known (v). Now I want to find a new distribution where i calculate the distance of all estimates from the known value (|v-e|). Since this new distribution is bounded by zero, the new distribution is no longer normal is seems. Does this make the new distribution equivalent to a Poisson distribution?


title: Checking for outliers in a glmer (lme4 package) with 3 random factors

text: I have a question relating to the checking for outliers and / or influential points in my dataset using a `glmer` model with 3 random variables. I'm investigating the detection rate (`SumDetections`) of receivers over increasing distance (`sc.c.distance`), and the effect of environmental influences on this (`depth`, `temperature` and `wind`) and how this differs between different transmitters used, controlling for random effects of `receiver ID`, `replicate` and `area`. I found that the `influence.ME` might be of help, so I checked it out. In the manual I read that this package is only able to delete levels of 1 single grouping factor or 1 data point per time over the whole data set. Unless I read the package info incorrectly, this package cannot do what I'm looking for. I'm looking for a way to check for outliers nested within 4 grouping layers. How my data is organized is as follows: First my data discerns between `Areas`. Within areas, multiple `replicates` were done. Each replicate consisted of 5 `distances` at which the detection rate was tested. For each distance, 20 `receivers` were tested. My model looks like this:               m <- lmer(SumDetections ~ tm + sc.c.distance + tm:sc.c.distance + c.tm.depth +                   c.receiver.depth + c.temp + c.wind + (1|replicate) + (1|SUR.ID) + (1|Area),                   data = df3, family = poisson)      My questions are:   * Is it possible to check for outliers of which the data is nested within 3 layers with use of influence.me?    * If so, how should I specify the command to get what I'm looking for, and how should I interpret the returned data by the Cook's distance or dfbetas?    * If not, is there another package that allows me to check for outliers?


title: What do you call a parameter that is estimated from historical values?

text: There are several methods to estimate parameters in a model (MLE, MAP, GMM). Does the process of estimating a parameter from historical data have a name?


title: How to calculate p values in logistic regression with gradient descent algorithm

text: In logistic regression, the gradient descent algorithm for calculating coefficients can be described this way: Until convergence, do $$ \beta := \beta + \alpha \frac{\partial L}{\partial \beta} $$ where L is the log likelihood function. Coefficients can be got quite easily here, but how can I get the p values with this algorithm?


title: regression vs correlation: explanatory power

text: I have a dataset where each instance is yes/no for a set of features, and has a score. Let's pretend they're sundaes. Thus [1,0,1]->6.5 Might mean yes chocolate, no vanilla, yes strawberry with a rating of 6.5 (I actually have 50 features and ~10,000 instances) I'm interested in seeing which features tend to produce higher ratings (like with sundaes, an instance can have any number of the features true, but typically between 2 and 5 are true). Currently I've tried (1) finding the correlation between score and each feature, and (2) performing linear regression and looking at the weights of the resulting model. What I don't understand is which is better. The first method gives me a certainty metric (I.e. p value) and the second method produces a generative model from which I can make predictions. But is one better than the other for saying which features are better choices? Is there a better option than either?


title: Cluster analysis on related factors

text: I am analyzing a public data set of information security incident data and trying to find "clusters" of related factors. Specifically, each incident is analyzed using VERIS for the actor's variety (e.g. actor.external.variety = hacktivist or actor.internal.variety = sysadmin) and motive (e.g. actor.motive = grudge) as well as the action's variety (e.g. action.hacking.variety = DoS or action.social.variety = phishing). The possible enumerations are pre-determined. In terms of size, I have a few thousand incidents and a few hundred factors / features, so it's large enough to make manual analysis difficult but small enough that runtime rarely becomes an issue. I want to find "clusters" that represent highly-correlated sets of factors, or incident patterns. As an example, we know that ideology motives and hacktivists and denial-of-service attacks frequently occur together. My current approach has been to create a sort of boolean matrix, where the rows are the individual incidents and the columns are the factors of interest. (This means I'm intentionally excluding things like asset data for now because those aren't the focus of my research.) Then, for each factor, I add up all the times it co-occurs with each other factor and build a weighted graph out of that to look for cliques and similar. Although it's been somewhat productive, I don't think this is the "right" approach. I also considered PCA for dimensionality reduction followed by kNN or similar for clustering. So my question is **what statistical / ML techniques best lend themselves to finding incident patterns in this type of data set?**


title: Where do I declare prior parameters in Stan?

text: I am fitting several multi-level hierarchical regression models using `rstan`. Each of these models has a unit information prior over the top-level distribution. The unit priors' parameters were originally store in an R `list`. My questions are:   1. Do I pass the unit information prior parameters for a particular model from the `list` where the rest of the data is stored?   2. Where do I declare these prior parameters? In the parameters block? Thanks.


title: Why does pre-rotation affect PCA?

text: So I noticed the following phenomena: I first generate a three-dimensional (x,y,z) multivariate normal distribution with no covariance between the variables. I then generate a second distribution by rotating the first distribution 90 degrees around the x-axis. I then perform PCA on each distribution to find the rotation matrix that will align them to the principal moments, and rotate each distribution accordingly. I thought that the original distribution and the rotated distribution should be identical after aligning with their respective principal moments, but I find that they continue to be perpendicular in the x dimension. Is there a reason for this? I'm not sure I can understand how pre-rotation will affect the results of PCA. To do this, I am using princomp in R. Thanks!


title: Q-Value Less than P-Value

text: I have a list of p-values from a series of regression equations. Many of these P-values are small. To compute false discovery rate I have been using the q-value package in R. I notice, based on this analysis, that P-values of 0.05 correspond to q-values that are smaller than the p values (0.03). Would I then say that my false discovery rate is smaller than my P value? Is this possible?


title: Best statistical test to show significance of two data coming from the same population?

text: Hope i'm posting in the right thread. I have a sample of genes where I calculate a number of features in the start of each gene, and a number of features in the end of each gene. I plotted the results in 2 curves, one is the number of features in the start of each gene, and the other one is the number of features in the end of each gene, according to the position in the gene. Now i'd like to know if the results I got are significant, knowing that the 2 curves come from the same population of genes (difference is the number of features is not calculated at the same location of the genes). What would be the best statistical test to see if my results are significant? Thanks in advance. PS: I found this about standard errors. Do you think in my case, I could just use standard errors and if they don't overlap between the 2 curves, i can say my data are significant? http://www.biochemia-medica.com/content/standard- error-meaning-and-interpretation EDIT: What i'd like to do is to show that i didn't get those curves by chance. The y axis represents the number of features and the x axis represents the position in the gene. I used the same sample of genes for the 2 measures plotted in the 2 different curves. Is there any test I can do to show that my results are good? And here the corresponding curves of my results with the standard errors for each one of them: ![enter image description here](http://i.stack.imgur.com/9QCdI.png)


title: Problem with lqmm. Very long time to run and issues with error handling

text: I am trying to use `lqmm` to fit some models and it appears to be working ok with my data (N 2000 w/ clusters=1200). I am using a call like this:               lqmm (y~x + covariates , random=~1,group=subjectID, iota=c(0.25, 0.50, 0.75),...)      However, I am getting these optimization errors every time:               Warning message: In errorHandling(OPTIMIZATION$low_loop, "low", control$LP_max_iter,     : Lower loop did not converge in: lqmm.     Try increasing max number of iterations (500) or tolerance (1e-05)       I have increased the number of max iterations to 500 and tolerance to 1e-05 as instructed but no dice. I still get the warning! Also, when attempt to obtain summary(object), it takes very very long time (3090 seconds for one that worked). Has anybody encountered these problems before and would like to share some ideas about how to go about this?


title: Case Control Study

text: Suppose we have identified $20$ cases from a population of $1000$ people. $15$ of these cases are male. We want to use one-to-one matching to choose the control group. We would match on gender. So would we just randomly choose $15$ males from the population that are not cases and then arbitrarily match them to the $15$ cases? There are $\binom{15}{2}$ matchings possible.


title: Model assumptions not met but model has predictive capabilities

text: This is a general question, has probably been asked before (I searched and did not see similar), may have a simple answer, may be ridiculous, and is just to satisfy my own curiosity. Say we build a linear model (or any model for that matter) on our training data, and then generate predictions from our test data and find that the model yields decent predictions as measured by say r^2. We then test the assumptions of our model, find that it fails to meet necessary assumptions, and are either unable to manipulate our data meet the assumptions or able to make our data to meet assumptions, but the new model built performs poorly compared to the invalid model. Any reason (besides that of a theoretical standpoint) not to use the better- performing invalid model, so long as it keeps generating accurate predictions?


title: jtest vs coxtest

text: In library "lmtest" The jtest and coxtest functions are used to see if two non-nested models are equivalent. According to the help function:               jtest performs the Davidson-MacKinnon J test for comparing non-nested models.           coxtest performs the Cox test for comparing two non-nested models.       Under what conditions should one be used over the other?


title: How to create a single rank / score number from a multiple set of variables that define the quality of something

text: I have several variables that describe the quality of a trading relationship between a client and financial institutions (FIs) the client deals with. The variables are populated with observed data over a period of time. Lets call the variables A,B,C,D. The client is able to use this data to filter on FIs based on either A,B,C or D right now. What I am trying to do is create a single rank/score which aggregates a chosen number of variables and a defined weighting to each variable based on what the client cares about. So for example, client ACME might create a rank based on xA + yB where x and y are weights applied to the two variables it cares about with other variables are ignored. I need this rank to be simple and easy to explain to clients who are not stats experts. I welcome advice and suggestions on how to approach this problem.


title: How would I factor study design into relative risk?

text: I'm comparing various conditions versus some disease for a large stratified/clustered dataset which purports to account for the entire population after weighting -- hence the use of RR instead of OR. An epidemiologist colleague suggested that I include relative risk in my analysis, but I feel that the version of the confidence interval available on wikipedia isn't sufficient, since it doesn't take sample design into account. Here's my SAS syntax, as recommended by the sample designers:               proc surveyfreq data=data.data_set_name;          tables disease*var1 disease*var2 disease*var1*var2;          weight WEIGHT_VAR; cluster CLUSTER_VAR;          strata STRATIFYING_VAR;          strata YEAR;      run;      My SAS output has the following columns I feel should be applicable: `Weighted Frequency` (the numbers I plug into the table), and `Std Dev of Wgt Freq` (actually standard error; SAS is weird). My intuition is that $$95\%CI=\exp(m \pm 1.96s) = \exp(\ln(RR)\pm1.96s)=\exp(\ln(\frac{\frac{a}{b}}{\frac{c}{d}}))\pm1.96\sqrt{\frac{b}{a(a+b)}+\frac{d}{c(c+d)}})$$ where $s$ is the square root of the variance, might be adjustable for the standard error in some way, maybe by applying the standard error to each value and propagating to the CI somehow. Edit: edited to include information mentioned in comments.


title: Glm with Poisson in R with no p-value?

text: I've been doing glm and thought I had it cracked but looking at some sites for interpreting results I realised R had displayed my results oddly.               > glm(Number.of.Fish~Julian.Date*Size.Class,na.action=na.omit,poisson)     Call:  glm(formula = Number.of.Fish ~ Julian.Date * Size.Class, family = poisson,          na.action = na.omit)     Coefficients:                       (Intercept)                    Julian.Date                          -0.0315824                      0.0003701                    Size.Classmedium              Size.Classmedium                            0.1247642                     -0.0268892                     Size.Classsmall   Julian.Date:Size.Classmedium                           3.6317021                     -0.0005089       Julian.Date:Size.Classmedium     Julian.Date:Size.Classsmall                                  NA                     -0.0156465       Degrees of Freedom: 355 Total (i.e. Null);  349 Residual       (298 observations deleted due to missingness)     Null Deviance:      2196      Residual Deviance: 1368         AIC: 2213      I've tried changing the code to family="poisson" but no change and I' afraid after that I'm stumped.


title: How do I test for homogeneity of regression in MANCOVA?

text: I am using a one-way independent MANCOVA with 4 dependent variables and a single covariate. My predictor is bivariate. I'm currently trying to test for heterogeneity of regression in SPSS by assessing the interaction between the covariate and predictor. The question is whether I should test the interaction for each dependent variable at a time or for the MANCOVA model as a whole? Preliminarily, the interaction in the MANCOVA model is non-significant, while two of of the interactions in my multiple ANCOVA comparisons are significant. Currently this leaves me at an impasse. I believe the the multiple ANCOVA interactions will inflate the familywise error of finding significant heterogeneity – is this assumption correct? Theoretically, I have reason to infer homogeneity of regression.


title: How to replace missing values with unsupervised random forest?

text: from http://www.stat.berkeley.edu/~breiman/RandomForests/cc_home.htm#unsup : > Formulating it as a two class problem has a number of payoffs. Missing > values can be replaced effectively. Outliers can be found. Variable > importance can be measured. Scaling can be performed (in this case, if the > original data had labels, the unsupervised scaling often retains the > structure of the original scaling). But the most important payoff is the > possibility of clustering. I tried searching, but couldn't find any additional information.


title: What does "(scores mapped into 1:m using rounded scores)" from the 1-sample permutation test mean?

text: I'm trying to do an exact paired permutation test. What does this "mapping" mean?               > library(exactRankTests)     > perm.test(rnorm(10), rnorm(10), paired=TRUE, exact=TRUE)              1-sample Permutation Test (scores mapped into 1:m using rounded scores)          data:  rnorm(10) and rnorm(10)     T = 6, p-value = 0.1836     alternative hypothesis: true mu is not equal to 0


title: Importance of original variables after PCA

text: I have a data frame in R containing four variables. I ran a PCA which showed that first two principal components mostly explain much of the variance in my data. I know that principal components are linear combinations of the original variables. Now how do I find out which combination(s) of my original variables give rise to these PCs? Thanks


title: Population Variance and Sample Variance

text: How is it that lowering the number of frequence,$n$ by $1$ in the for formula for Population and Sample Variance account for the discrepency of using sample rather than population. I mean how is dividing the whole expression by $n-1$ in sample variance better than dividing by $n$? I am referring to the formula $$s^2=\frac{\sum_{i=1}^{n}({x_i- x_{avg}})^2}{n-1}$$ more accurate than $$s^2=\frac{\sum_{i=1}^{n}({x_i- x_{avg}})^2}{n}$$ Thank You.


title: Finding a subset of the data in which two variables are independent

text: I have a dataset of ~400 subjects. For each subject, I have about 18 variables (V1...V18), all of those continuous and measured at the same timepoint. A paper I'm working on would benefit from showing that one of those variables (say V1) has a significant correlation with one of three other variables (say V2, V3, or V4) – all of these are measuring different aspects of the same illness. And I have found some significant correlations between the pairs (V1–V2, V1–V3, V1–V4). However, I have another variable, say V5, that is more strongly correlated to all of them (V1–V4). When I correct (linearly) for V5, the correlations between V1–(V2/V3/V4) is no longer significant/there. So at this point, any critic could say that the V1 is correlated to V5, V5 is correlated to V2/3/4, and that I can't make the point I want to make. However, I have reason to believe that subgroups of the subjects (more susceptible to the illness), will show strong correlation between V1 and V2–4 that is independent of V5. _Is there any way to find these subgroups?_ I would like to find some partition(s), A < Vx < B, in which there is a significant correlation between V1 and V2/V3/V4, that I can show is independent of V5. Alternatively, if there's a way to quantify how much of the correlation between V1 and V2–4 is due to correlation between V1–V5 and how much is independent, that would be a good solution as well. I have been reading a lot about two-way ANOVA, thinking about doing it between V1 as rows, V5 as columns, and one of V2–4 as outcome. But even if I do this (which would take significant effort, since I couldn't find a simple Python implementation), I'm not sure if it would advance me towards my goals. Any help is appreciated.


title: Estimating how precise a linear regression is

text: I just started learning regression, so my background is not very strong. I have to solve a certain problem, which doesn't seem too hard, however I have some trouble in understanding it. Suppose we define               x<-seq(1,50,1), y<-3*x+2+6*rnorm(length(x))      Hence my response $Y$ is exactly of this type: $Y=X\beta+\epsilon$, where $\epsilon$ is the error and normally distributed (here with mean 0 and variance 1). Now we should do a simple linear regression $n$ times and saving the $n$ slopes of the regression.               reg <-lm(y~x)      The true value of the slope should be $3$, as assumed. I calculated the $n$ slopes and saved them in a vector $s$. The next task is to draw a histogram of the $n$ estimated slopes and add the normal density of the theoretically true distribution of the slopes. Here my problem starts: it is clear that the slopes are normally distributed. Hence I have to calculate the mean and variance. Looking at the model: $y_i=\beta_1+x_{i,2}\beta_2 +\epsilon_i$, where $\epsilon_i\sim\mathcal{N}(0,6^2)$. The hint is to use `solve` for the inverse of a matrix. I do not see how this should help? Is it meant, that I should calculate this: $$X^{-1}(Y-\epsilon)=\beta$$ How do I get from there the mean and variance? And what is the difference to the mean and empirical standard deviation of the estimated slopes? Thanks in advance for your help.


title: Defining PSU in "Sampled with Replacement" Cluster Samples

text: I am trying to decide on an approach to estimate design effect for a multi stage cluster survey. The clusters were selected with probability proportional to size sampling WITH replacement. The primary sampling units (districts) are large enough that some were sampled multiple times. When defining the PSU variable should the unique districts within the sample be used, or a cluster variable representing the true number of clusters sampled? For example: A district gets sampled twice, should the PSU variable described, be one or two individual clusters? Thanks!


title: Which clustering algorithm shall I use?

text: I need some help My project aims to develop algorithms for spatial temporal analysis of Flickr, Twitter and Foursquare databases to detect any kind of significant changes, named as “Event” in real time. Event can be defined as any anomalous user activity, which happens at a time or within a specific period of time. For this, different clustering methods should be implemented and the best fit has to be selected. The detected events will be visualized for further exploration. This information will be integrated with some other VGI sources to provide a series of Volunteered Geographic services. Could you please suggest me which clustering algorithm is good for this project? And also please suggest me some books and study material...


title: Estimated marginal means and means from descriptive statistics - interpretation in ANCOVA?

text: I ran an ANCOVA which proved that the covariate is significant. If i understood right, the estimated marginal means show what are the means of the DV when we account for the covariate. Then we could compare them to the means from the descriptive statistics and see what is the difference under the influence of this covariate. In my case, the covariate is "brand commitment" and the DV is "purchase probability". So what I get is - when I have large value for the brand commitment, the mean purchase probability increases for all my experimental conditions and vice versa when I have low. If someone can confirm that my observations are right, I will appreciate it. Also, is there a specific way to report this i.e. to sound more "academical"? And is there something else important for interpretation to look at from the ANCOVA results, besides the f values, parameter estimates? Thanks.


title: Maximum likelihood hypothesis vs Expectation maximization

text: Maximum Likelihood is given by the formula $h_{ML}=arg\space max_{h \in H } \space\space p(D/h)$ I want to transform it in terms of mechanism involved in Expectation maximization. $h_{ML}=arg\space max_{h \in H } \space\space p(D/h)$ $\implies h_{ML}=arg\space max_{h \in H } \space\space \ln (p(D/h))$ $\implies h_{ML}=arg\space max_{h \in H } \space\space \ln (p(D/h))$ Now i divided the above in to two parts . $Q(h'/h) \leftarrow E_{h,X}[\ln (p(D/h))]$ $h_{ML}=arg\space max\space Q(h'/h)$ $Q$ is just an intermediate function. $E_{h,X}$ is an expectation function that expects unknown data w.r.t observed data $X$. Is it correct way of interpretation ?


title: How to analyse data when it comes from different conditions

text: In my experiment, I am varying 4 different parameters (`side`, `eccentricity`, `latency` and `duration`), and for each of these conditions I have 2 or 4 values (e.g., `eccentricity` can be `5 deg` or `10 deg`, `side` can be `right`, `left`, `bottom`, `top`, and so on). I want to compare the rates for `5 deg` vs. `10 deg` and for both of them, the data can be either from `top`, `bottom`, `right`, `left`. My concern is that if I just add rate across all sides, the `side` condition will affect my `eccentricity` condition and my error will be quite high. So if I cannot sum the data from the sides, how can I compare the means of `5 deg` vs. `10 deg`, while taking the sides into account?


title: Understanding the derivation of a ML-estimator

text: I have the following likelihood function: ![enter image description here](http://i.stack.imgur.com/1d5Tw.gif) I'm given this information about the $\Omega$ matrix ($\boldsymbol{1}$ is a $T \times 1$ vector of ones): ![enter image description here](http://i.stack.imgur.com/YLjFa.gif) I would like to be able to show that the likelihood function can be rewritten to this: ![enter image description here](http://i.stack.imgur.com/sKuPa.gif) I think there is a misprint in the article, hence the screen images from the article I attached are wrong, meaning that $\Omega^{-1}$ is supposed to be $$\Omega^{-1}=\frac{1}{\sigma^2_{\varepsilon}+T\sigma^2_{c}}\boldsymbol{1}(\boldsymbol{1}'\boldsymbol{1})^{-1}\boldsymbol{1}'+\frac{1}{\sigma^2_{\varepsilon}}(I_t-\boldsymbol{1}(\boldsymbol{1}'\boldsymbol{1})^{-1}\boldsymbol{1}')$$ and not as stated above. So far I've inserted the expression for $\Omega^{-1}$ in the top log- likelihood function. However, especially the last expression of the rewritten log-likelihood function bothers me, i.e. this part: \begin{equation*} \begin{split} ...& &+\frac{1}{2}\frac{1}{\sigma^2_{\varepsilon}}\frac{1}{N}\sum_{i=1}^{N} \left(Y_i-\boldsymbol{1}\bar{Y}_i-\left(X_i-\boldsymbol{1}\bar{X}_i\right)\beta\right)'\left(Y_i-\boldsymbol{1}\bar{Y}_i-\left(X_i-\boldsymbol{1}\bar{X}_i\right)\beta\right). \end{split} \end{equation*} Does anyone see the solution? I also have this lemma, but not sure if it is to be used in the rewritting of the likelihood function: ![enter image description here](http://i.stack.imgur.com/wlhbU.gif)


title: Choosing the right regression and validating it

text: I have two questions:   1. _Which regression should I choose and why?_ Based on the R squared value, exponential regression seems to be a better fit. But I am not really sure, if I should just go with it because it has a higher R squared value? ![enter image description here](http://i.stack.imgur.com/2XIPi.png)   2. _How do I validate it?_ Lets say I decided to go with the linear regression (despite its smaller R squared value.) The only process I am following is to manually calculate the residuals and fit a histogram to check if they are normally distributed. And apparently they are OK in this case. However, my concern is that some values (red circled ones) don't seem to be fitting well since their error is too high. But still, I suspect that they do not undermine R squared value, because are symmetric to each other. So it will probably add a bias? So what is to best way to validate my choice of regression against such cases (besides high R squared value and normality check.) And are there other biases that I should look for?


title: Modelling uncertainty with probability distributions

text: Suppose we’re building a system with certain qualities that we’re interested in, e.g. response time, battery usage, etc. Each of the system’s qualities depends on our decision about components we use to build the system. We only know how the **decisions** impact the **qualities** of our system with uncertainty, e.g. deciding to use a GPS module will most likely add an extra 10µJ to our battery usage, but could add anything between 9-14µJ. Is it common to model this uncertainty using probability distributions? Is beta distribution a good choice here? Or would triangular or uniform distribution be more sensible in this case? In case of beta distribution, how do I work out the α and β parameters knowing the most likely and limiting estimates only? Note: my knowledge in statistics is rather limited.


title: Merge categories or drop one?

text: On seeing the following table, I concluded there is a lot of sparseness in the 'very unsafe' category. My first reaction was to merge it with the 'unsafe' category, leaving me with 3 categories. But another option is to just delete the category. What are the considerations to make before choosing one option? Also, I am wondering if it is ok to keep 2 degrees of safe, and only 1 degree of unsafe. It would, to me, not make much sense when discussing the odds. ![enter image description here](http://i.stack.imgur.com/vSNvO.png) edit in response to Nick Cox: I assume I should remove the category since the conclusions seem out of line with theoretical expectations. For example, it is assumed that those living in more urban areas experience greater levels of fear. And indeed, the odds of feeling very safe (vs safe) are: countryside > small city > big city. The odds of feeling safe (vs unsafe) have the same order. But then for unsafe (vs very unsafe) it becomes countryside > big city > small city, which is out of line with theory and the former trend. I could say that, although urbanisation brings with greater feelings of unsafety, these feeling take on less severe degrees when initiated urbanisation proceeds to higher levels. But of course, I don't want to report such conclusions when these are due to data sparseness rather than social realities.


title: How to do 4-parametric regression for ELISA data in R

text: I am a biology student. We do many Enzyme Linked Immunosorbent Assay (ELISA) experiments and Bradford detection. A 4-parametric logistic regression (reference) is often used for regression these data following this function: $$ F(x) = \left(\frac{A-D}{1+(x/C)^B}\right) + D $$ How can I do this in `R`? I want to get the $A$, $B$, $C$ and $D$ values and plot the curve. PS. If I have some data, how can I use the calculated function $F(x)$ to get the value? I mean how do I go from "data -> F(x) -> value"?


title: transformed data- effect sizes

text: I have been working on some data in collaboration with a statistician. As the data was not normally distributed a log-transformation was performed. However, it was explained to me that we cannot report effect sizes of our findings in the original units because of this transformation. We therefore resorted to using percentages. I need to explain the reason for this in simple terms to collaborators (who are in the clinical field, not statisticians). Does anyone have suggestions of how to do this?


title: How to assess correlation when each variable is measured by independent replicates?

text: I frequently measure multiple variables, in multiple replicates, at many sites. For instance, I might measure bacterial abundance and bacterial growth rates, each in 3 replicates, at many sites. Each replicate is independently sampled _and_ each variable is measured in a different sample (i.e., I can't measure both bacterial abundance and growth rate in the same sample). I'd like to test for correlation between those variables. The problem is, since the variables are measured independently, the variables are not paired. Replicate 1 of variable A is not related to replicate 1 of variable B, any more than it is to replicate 2 of variable B. I could test for correlation among the average of replicates at each site - but that seems blunt, since you lose information about variation among replicates for each parameter. I can imagine some kind of resampling approach, where I randomly select one replicate for each variable at each site. Is there a better way?


title: Identity of moment-generating functions

text: Are there any non-identical distributions which happen to have the same moment-generating function?


title: Comparing multiple regression to iterated linear regressions

text: Will running multiple regression including 2 independent variables give the same output for each independent variable if we run linear regression using the same independent variables separately/individually?


title: Plot result of Nemenyi's test

text: Could you please tell me which R function/library (not necessarily R) to use to plot Nemenyi's test result with the critical differences? Below there is an example which was retrieved from this* paper (page 9).      ![enter image description here](http://i.stack.imgur.com/SPp8E.png) ***Reference** : BROWN, I. An experimental comparison of classification techniques for imbalanced credit scoring data sets using SAS® Enterprise Miner™. SAS Global Forum 2012. Paper 129-2012. University of Southampton, Southampton, UK.


title: Interpreting R results, are the data multivariate normal?

text: I ran "mvn" using the "mclust" package in R using the following codes:               mvn("EEE", data[,18:22], prior = NULL, warn= NULL)      I am having trouble figuring out how to tell if the results are indicating "multivariate normal fit". The results are below:               $n     [1] 707          $d     [1] 5          $G     [1] 1          $parameters         $parameters$pro     [1] 1          $parameters$mean             [,1]     V18 2.993399     V19 3.775578     V20 3.790193     V21 2.652051     V22 3.235266          $parameters$variance     $parameters$variance$modelName     [1] "XXX"          $parameters$variance$d     [1] 5          $parameters$variance$G     [1] 1          $parameters$variance$Sigma                 V18          V19          V20          V21         V22     V18  1.11483935 -0.046114371  0.021088813  0.094512987 -0.03160755     V19 -0.04611437  0.604042865 -0.001037869  0.001748973  0.13845048     V20  0.02108881 -0.001037869  0.511536706  0.112759387 -0.17223313     V21  0.09451299  0.001748973  0.112759387  0.542454933 -0.13187494     V22 -0.03160755  0.138450479 -0.172233126 -0.131874937  0.81860855          $parameters$variance$cholSigma             [,1]        [,2]          [,3]         [,4]        [,5]     [1,] 1.05586 -0.04367472  0.0199731239  0.089512842 -0.02993537     [2,] 0.00000  0.77597383 -0.0002133427  0.007292026  0.17673671     [3,] 0.00000  0.00000000 -0.7149389730 -0.155220369  0.24001701     [4,] 0.00000  0.00000000  0.0000000000  0.714349947 -0.13050827     [5,] 0.00000  0.00000000  0.0000000000  0.000000000  0.84370373          $parameters$variance$sigma     , , 1                      V18          V19          V20          V21         V22     V18  1.11483935 -0.046114371  0.021088813  0.094512987 -0.03160755     V19 -0.04611437  0.604042865 -0.001037869  0.001748973  0.13845048     V20  0.02108881 -0.001037869  0.511536706  0.112759387 -0.17223313     V21  0.09451299  0.001748973  0.112759387  0.542454933 -0.13187494     V22 -0.03160755  0.138450479 -0.172233126 -0.131874937  0.81860855                         $loglik     [1] -4279.837          attr(,"returnCode")     [1] 0


title: Estimating variability in a response variable not accounted for by measured predictors using linear regression and R

text: ### Background I have a data set of patients who were operated on at two different hospitals, A and B. Lymph nodes were removed from each patient during the operation and counted, this is saved as `LN_reviewed` for each patient. I want to know how much variability there is in the lymph node number that is **not** accounted for by gender, the year of the operation, or the age of the patient when operated on. My assumption (hypothesis) is that any additional variability is likely due to the pathologist at the institution (this was not actually measured in my study). ### My question What is the best way to go about estimating the variability in lymph node number that is not accounted for by gender (a factor), year (a continuous variable), or age (a continuous variable)? ### My initial attempt at answering this question I built a linear regression model using the number of lymph nodes as the response variable and the gender, operation year, and operation age as predictors. I am not sure how to interpret the results to answer my specific question. Should I be looking at the R squared? If so, is there a way to get a confidence interval for it? Thanks to anyone who can help. If you think I am going about this the wrong way, please let me know.               Call:     lm(formula = LN_reviewed ~ Gender + Operation__year + Operation__age, data = sample_data)          Residuals:         Min      1Q  Median      3Q     Max      -49.436 -15.280  -0.495  13.450  61.564           Coefficients:                       Estimate Std. Error t value Pr(>|t|)         (Intercept)     -1.190e+04  1.459e+03  -8.159 9.95e-14 ***     GenderMALE      -5.542e+00  4.685e+00  -1.183    0.239         Operation__year  5.980e+00  7.296e-01   8.196 8.01e-14 ***     Operation__age  -2.524e-01  1.675e-01  -1.507    0.134         ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1          Residual standard error: 22.46 on 158 degrees of freedom     Multiple R-squared:  0.2999,    Adjusted R-squared:  0.2866      F-statistic: 22.56 on 3 and 158 DF,  p-value: 3.268e-12


title: Linearity between predictors and dependent variable in a linear model

text: I run the following linear model in R :               lm(formula = NA. ~ PC + I(1/SPCI), data = DSET)      The p-value for each predictor is significant, and it works fairly well with respect to most of the assumptions in linear regression, such as:   * Normal distribution of errors.   * High correlation between predicted values and estimated values.   * Homoscedasticity.   * Non-collinearity between the predictors: PC and (1/SPCI) are not correlated at all! But, digging more into the topic, there's an assumption that fails in my model, and it says: > Linearity of the relationship between dependent and independent variables. This kind of contradicts the non-collinearity assumption because, if NA. and PC are highly correlated and NA. and (1/SPCI) are too, then PC and (1/SPCI) are correlated, and this is violating the assumption of non-collinearity. I think that I misunderstood the assumption, or there's an explanation to this.


title: Guessing the correlation structure of a mixed model

text: If you would fit a small sample like with the standard framework for linear models. Would you suggest from the residuals any special correlation structure like compound symmetry or first order autoregression?               Y<-data.frame(response=c(10,19,27,28,9,13,25,29,4,10,20,18,5,6,12,17),                    treatment=factor(rep(1:4,4)),                    subject=factor(rep(1:4,each=4))                    )     fit<-lm(response~-1+treatment,Y)


title: opinion poll, t test or paired test or others

text: Which test should I use to test the difference of monthly opinion poll? If I want to know whether the poll result of April is significantly different from March, should I use independent t test or paired test? If I use independent t test, it seems that I am comparing two groups of people instead of change over time, however the samples be totally different between two time points. So, I am not sure which to use. Thx!!


title: Found a better measurement - what to do with the old/new measurement (both formative)?

text: In a survey, I used a measurement of uncertainty as a key dependent variable, and after implementing the measurement in 2 waves of my survey, I found a better measurement and this new measurement has better internal validity and fits better theoretically. I still have 2 waves of survey to be sent and I want to add the new measurement. How to proceed with the analysis? I will likely have 2 waves of data with the old measures, and 2 waves with both the old and the new measures. The waves of data are just batches of sampling, but it seems all 4 waves are needed to have a good enough sampling size. Thank you!!


title: How can I simulate data for comparing logit and probit model?

text: Can anybody tell me how can I simulate a data set for comparing logit and probit model? If I simulate a data from a standard multivariate normal distribution considering 5 variables then can I use one of them as continuous latent variable (for which I consider non-zero covariance in the covariance matrix)? If I can, then what should be the cut point (of latent variable $y^{*}$) for considering $Y=1$ and $Y=0$. And I need some reference where this short of works are conducted.


title: Logistic regression understanding variable nature

text: I have 2 categorically dependent variables(both binomial) in logistic regression which individually both give positive estimates against the response(binomial). However if modeled together one give positive and other give negative. I checked for the multicollinearity using VIF, which seems fine. What could be the possible reason for such behaviour.


title: Computation of parameters of truncated normal distribution

text: I have a numeric data sample with mean $\mu$ and standard deviation $\sigma$. I believe that this data has a normal distribution truncated to $[0,1]$.   1. Is there a reasonably simple formula for estimating the parameters $\mu', \sigma'$ of that truncated normal distribution? Is there an R command which does it?   2. I understand that I could find $\mu', \sigma'$ by applying a maximum-likelihood fitting algorithm directly to my data (through, say, fitdistr in R).   3. Is there some good understanding how much more precise approach (2) is over approach (1)? Say the sample size is from 10 to 100, $\mu=0.75,$ $\sigma=0.15$.


title: Comparing dependent correlations

text: I have a repeated measures data set with each subject ($N = 32$) being measured with two repeats of two tasks at two levels. My values of the correlation coefficient $r$ are:               1st repeat vs 2nd repeat: 0.6847     level 1 vs level 2: 0.6930     task A vs task B: 0.4081      The assumption is that the 1st and 2nd repeat of a task at a particular level should be noisy measurements of the same thing. I think, what I would like to know is if task A and task B are measuring different things and if a task at level 1 is measuring the same thing as a task at level 2. I think this might be similar to testing if the repeat correlation is different from the task/level correlations. I guess my first question is, is this true. My second question is, how do I compare correlation coefficients with a repeated measures design?


title: Permutation Analysis using unequal samples

text: I need to conduct a permutation analysis on unequal samples, one (X) with 579 data points, and another (Y) with 1289 data points. There is an obvious difference in the data. However, this difference may be due to the imbalance between sample sizes. I need to run a permutation analysis using all possible pairs (combinations) of these data points to see what the likelihood that these differences are real and not artifact due to the imbalance. The results need to give me p-values to interpret for significance testing. I am very new to R. Also, most everything I have found gives me paired-samples. Does anyone have code or suggestions to do this? Thanks.


title: Choosing between one and two clusters using Calinski-Harabasz criterion

text: Using R: `pamk.best <- pamk(data,k=2:12,criterion="ch")` returns an estimation of the best number of clusters using Calinski-Harabasz criterion. However, this criterion cannot be used to evaluate the quality of 'one' cluster. Therefore, how can I know if it is better to divide the data into two groups or to keep it unchanged?


title: Find the distribution parameters

text: On STAT101 your assessment is based on: Final Exam 50% Learn based on-line assessment 30% Assignments 20% Consider three random variables $X$, $Y$ and $Z$ which respectively represent the exam, on-line assessment total and assignment scores (out of 100%) of a randomly chosen student. Assume that $X$, $Y$ and $Z$ are independent (this is clearly not true, but the answers may be a reasonable approximation). Suppose that past experience suggests the following properties of these assessment items (each out of 100%): $E \left ( X \right ) = 61$, $sd(X) = 20$, $E(Y) = 72$, $sd(Y) = 22$ and $E(Z) =65$, $sd(Z) = 24$. a) Find the distribution parameters, $E(T)$ and $Var(T)$, for the total mark, $T$, where: $T = 0.5X+0.3Y+0.2Z$.


title: Pseudo-likelihood and RBMs

text: I need to train a restricted Boltzmann machine to model the joint probability of categorical variables. For this I adapted a Bernoulli RBM to have groups of softmax units in the visible layer. The original scikit-learn code I started with trains a Bernoulli RBM and computes the pseudo-likelihood value at each training iteration to measure progress. My question is: **Having changed the probabilistic interpretation of the visible units of the RBM (Bernoulli RBM to Softmax RBM), is it correct to continue using the same pseudo-likelihood computation as before? Or should it also be modified in some way?** Since the conditional independence between visible units despite a change in their probabilistic interpretation still holds, I'm inclined to think that it's ok to use the same pseudo-likelihood computation. However, I'm more of a programmer and my theoretic understanding of statistics is quite poor so I don't have a theoretical justification for this. I started reading the paper Besag, J. (1975) "Statistical Analysis of Non-Lattice Data." The Statistician, 24(3), 179–195 which seems to be a standard reference when it comes to pseudo-likelihood but I was also wondering if someone here could give me a simpler answer. Thanks!


title: Calculation the Expectation of an Inverse Wishart matrix

text: I have $\boldsymbol{A} = \boldsymbol{G}^H \boldsymbol{G}$ is a Wishart matrix, i.e, $\boldsymbol{G}^H \boldsymbol{G} \sim \mathcal{W}_K (M, \boldsymbol{\Lambda})$ with $\boldsymbol{\Lambda} = \mathrm{diag} (\rho_1, \cdots, \rho_K)$. ($\boldsymbol{G}$ is $M \times K$ matrix ). We know that $\boldsymbol{B} = (\boldsymbol{G}^H \boldsymbol{G})^{-1}$ is an Inverse Wishart matrix. So how can we calculate the expectation of $\boldsymbol{B}$, i.e, $\mathbb{E} \\{\boldsymbol{B}\\}$? Thanks!


title: How to assess the proportional hazards assumption for a continous variable

text: I am having a problem with checking the assumptions for a continuous variable in a proportional hazards model. If a variable were a factor with many levels, then I could use the logrank test or check whether the log(-log) transformations of survival curves are parallel. But what if a variable is continuous? Is that method still valid? Is Schoenfeld's test is a solution?


title: How can I create a testable hypotheses to compare benchmark results for the same tests run in 3 environments

text: I am planning on running the same benchmark tests (i.e. throughput, response time,..) in 3 environments: 1) a mobile device, 2) a virtual machine in VirtualBox, 3) a virtual machine on Amazon EC2. What could be the best way to define the hypothesis, Null or Alternate when considering significant differences or relations between the means obtained from the test runs? Any advice is welcome


Do the pdf and the pmf and the cdf contain the same information? For me the pdf gives the whole probability to a certain point(basically the area under the probability). The pmf give the probability of a certain point. The cdf give the probability under a certain point. So to me the pdf and cdf have the same information, but the pmf does not because it gives the probability for a point `x` on the distribution.


title: multivariate control chart

text: I am trying to build a multivariate control chart but I did not find any good explanation of how to do that. Starting with $x_0$ a part of the samples that are in control I compute $\mu_0$ and $\Sigma_0$ the mean and the covariance matrix. Now with new data $x_n$ I compute $d= (\bar{x}_n-\mu_0)\Sigma^{-1}(\bar{x}_n-\mu_0)$ . After this step I do not really know how to go further. The best answer will be a theoretical explanation and an R or python code.


title: Panel data and selection models issue

text: Good day,   I'm working with a panel dataset, I've used many models, homogeneous (fixed effect, pooled ols and Driscoll and Kraay) heterogeneous (swamy random coefficients) and would like to do a post-estimation to select the model that best fit my regression. Is there any method, command that may allow me to do this? Any hint will be highly appreciated.   Ama


title: Omitting a moderator

text: I'm wondering what the effect is if I don't include moderators in my model? Is this the same or different from an omitted variable bias? I am having a hard time grasping this conceptually. More information: I did a GEE analysis for a binary outcome for a publication. The independent variables were the study group and several demographic variables along with a few other participant characteristics. The CI of the odds ratios were quite wide, but I think it is due to the small sample size. One reviewer indicated we should include moderators in the model to decrease the variability of the estimates. I have never heard of this, so I was wondering if this was correct.


title: Can the formula developed by Olkin und Siotani (1964) also be used to compare correlations from Spearman's $\rho$ and Cramér's $V$

text: The followgin formula is normally used for comparing two Pearsons $r$ correlations: ![Formula from Olkin und Siotani \(1964\)](http://i.stack.imgur.com/hxL5j.png) Can this formula also be used to compare correlations from Spearman's $\rho$ and Cramér's $V$?


title: Which feature selection method to use for classification problem

text: I have to do some feature selection for a classification problem with numeric features. I am not sure which feature selection method to use. Chisquared test or Spearmann's rank correlation coefficient, which one would be the best. Any suggestion guys?


title: Paired or unpaired Wilcoxon test?

text: I counted the number of geese on an intertidal mudflat on 100+ days over the winter. I made two counts on each of these days: one at low tide and one at high tide. I want to know if the number of geese present differs at high and low tide. As the data are very positively skewed, using a Wilcoxon test is appropriate. However, should I use rank sum test (unpaired test) or signed rank test (paired test)?


title: Help with lmer in the lme4 package

text: I am new to using linear mixed models and would greatly appreciate any help I can get. I have an equation of the form $ y = X\beta + Zu + \epsilon$ where $u$ is a random effect whose covariance matrix I want to be a multiple of the identity matrix. My matrix $Z$ has about 1000 columns. Trouble is I have only one group, i.e. not multiple observations per person, and I think `lmer` does not like this. For example in the simple case, if I put in something like `lmer` $({\rm height} \sim 1 + (0+Z_1|ID) + (0+Z_2|ID))$ where $Z_1$ and $Z_2$ are 2 columns of $Z$ and $ID$ is the **unique ID for each person** `lmer` throws up an error saying I need the number of groups to be less than the number of observations (a similar problem occurs if I set everyone belonging to the same group). I had partial success with `lme` (in the `name` package) where I could provide the covariance structure of $u$ as in `pdDiag`, and the package allows one group. Trouble there was that I ran into memory issues, and the time was too long. Since `lmer` uses a lot of sparse matrix methods, I am guessing my flexibility in the **covariance structure** will be a lot lesser. That being said, my restriction seems like (probably?) one of the more obvious ones which would have been implemented in these models as is.


title: Number of observations in plm package, n, T, N

text: I am new to R. I am using the `plm` package to run pooled OLS and Fixed Effect. My data set has 130 countries over 33 years with lots of missing values. My model has 9 independent variables. When I run my models I get the following:               Unbalanced Panel: n=34, T=1-6, N=61      Can any body explain this? Does it mean that 61 countries have observations for all 33 years across the 9 independent variables? I saw the same question asked here: Interpreting plm output in R - number of observations used with very unbalanced panel Please help if you can. Thanks.


title: Cluster Analysis Sample Size for Diss Proposal

text: My psychology dissertation will be a cluster analysis with one grouping variable. How do I know how many participants I need? I imagine there should not be more than 5 clusters. 79 items make up 9 attributes that the participants will be clustered on. I realize that power analysis is moot since cluster analysis is not inferential. I have seen that for factor analysis it is 10 per item but cannot seem to find a guideline for cluster analysis. Thanks for any help or insight you can offer.


title: equivalent of PCA explained variance ratio for SVD ?

text: i am wondering if there is an equivalent of PCA explained variance ratio for SVD. What are the measures I can get to monitor the number of columns I keep after the SVD ? Are any of these metrics implemented in sklearn/python or R ?


title: Approximation for linear combination of truncated beta variables

text: Let $x_1, x_2$ be random variables that follow truncated beta distributions: $$f_{i}\left(x_{i}\right)=\frac{x_{i}^{\alpha_{i}-1}(1-x_{i})^{\beta_{i}-1}}{\mathrm{B}(a_{i},b_{i};\alpha_{i},\beta_{i})}$$ for $a_i\le x_i \le b_i$, and $f_i (x_i) = 0$ otherwise. Here $\mathrm{B}(a,b;\alpha,\beta)=\int_a^b t^{\alpha-1}(1-t)^{\beta-1}\mathrm{d}t$ is the incomplete beta function, and $0\le a_i \le b_i \le 1$ for $i=1,2$. Define the random variable $z = k_1 x_1 + k_2 x_2$, where $k_i$ are real constants. **Is there an approximate formula for the distribution of $z$?** A related problem is the linear combination of full beta distributions (that is, the particular case when $a_i =0$ and $b_i = 1$ for all $i$). A paper where this is treated is [1], where two approximations are studied. The first approximation is very simple, and suggests that $z/\sum_i k_i$ is approximately a Beta random variable with mean and variance equal to the exact mean and variance of $z/\sum_i k_i$. This is a very convenient approximation because the mean and variance of the linear combination are very easy to calculate, and the parameters of the Beta distribution ($\alpha$, $\beta$) are also easy to get from the mean and variance. I haven't found any more publications on this problem. The same idea seems to work well with linear combinations of Gamma distributions. That is, the linear combination of Gamma random variables distributes approximately as another Gamma variable with the mean and variance of the linear combination [2]. **Can we approximate a linear combination of _truncated_ Beta random variables (such as $z$ above) with a random variable that is proportional to a Beta random variable with the mean and variance of the linear combination? How good is this approximation?** [1] Jóhannesson, B., & Giri, N. (1995). On approximations involving the beta distribution. Communications in Statistics-Simulation and Computation, 24(2), 489–503. [2] Stewart, T., Strijbosch, L., Moors, H., & Batenburg, P. van. (2007). A simple approximation to the convolution of Gamma distributions.


title: How can I test for within-subject "consistency" across repeated measures?

text: I have a dataset with a clear question I want to answer, but I'm unsure of the best statistical approach to answer my question. Here's an example that matches the structure of my data: Let's say we have a group of ~200 physicians for whom we measure 10 comparable quality measures (e.g. adherence to guidelines), unit of measurement is a percent adherence to the measure. The question I want to ask is whether or not physicians tend to have consistently high or low adherence to all measures, or whether or not physicians have uneven adherence and tend to perform well on certain measures but not others. Ideally, the test or method I use can give me meaningful information for each individual. I think that the 10 separate quality measures can be conceptualized as a "repeated measure" because they are clinically very similar and I except some within-subject correlation. I pondered repeated-measures ANOVA to test the null hypothesis that there is no significant difference across the quality measures; but I worry that it will be too sensitive and almost surely reject the null hypothesis and I won't learn very much, especially at the individual physician level. I'd like to assign each physician to a percentile for the measures and perform some non-parametric test, but I'm getting out of my depth and was wondering if the CrossValidated community can help point me in the right direction. Suggestions for R packages also appreciated. Thanks!


title: Should the ridge parameter qualitatively change the weights?

text: I'm running a Ridge regression with about 1200 parameters (and about 30000 datapoints). I noticed that for some values of ridge, the weights look qualitatively different beyond a certain point. In particular, higher values for the ridge parameter will yield weights that clearly have some kind of structure to them, whereas at some point they look almost like white noise. In this case, I'd expect the weights to have the former kind of structure due to the features I'm using, but cross-validation often chooses a parameter for which the weights look totally noisy. It seems to me that the ridge parameter should only change the relative distribution of the weights, but not their qualitative differences. Am I wrong here?


title: Measure of association for 2x3 contingency table

text: I have a 2x3 contingency table - the row variable is a factor, the column variable is an ordered factor (ordinal level). I'd like to apply either symmetrical or asymmetrical association technique. What do you recommend me to do? Which technique do you find the most appropriate?


title: Implementing the 'kernel trick' for a support vector machine in R

text: I've heard a bit about the 'kernel trick' for support vector machines, and I was wondering:   1. How do you identify problems that might benefit from the kernel trick?   2. How to implement it in R? Thank you


title: Removing empty columns from a dataset

text: I am taking a class in data mining and I am working on a term project using the BRFSS dataset. I have a huge dataset with 405 columns and 12,000 rows. There are many columns which are completely empty. I was trying to remove empty columns using SAS, R or Excel but it doesn't work. Could you suggest a method to remove the empty columns or any tutorial that will help me cleaning up the data? There are a lot of missing cells too. I am using KNIME to train my data and it doesn't work if there are missing values. How can I handle the missing values?


title: Pearson correlation of data sets with possibly zero standard deviation?

text: I am having a problem computing the pearson correlation coefficient of data sets with possibly zero standard deviation (i.e. all data has the same value). Suppose that I have the following two data sets:               float x[] = {2, 2, 2, 3, 2};     float y[] = {2, 2, 2, 2, 2};      The correlation coefficient "r", would be computed using the following equation:               float r = covariance(x, y) / (std_dev(x) * std_dev(y));      However, because all data in data set "y" has the same value, the standard deviation std_dev(y) would be zero and "r" would be undefined. Is there any solution for this problem? Or should I use other methods to measure data relationship in this case?


title: Standardizing a non normal dataset

text: I will be begin by saying that I am not very familiar with statistics. I am calculating a drought index and I have a set of values. The next step in the procedure is to standardize the values to allow for comparing across regions and time scales. The paper that establishes this particular drought index fits their dataset to the log-logistic distribution and uses the probability distribution function to standardize the values. My dataset does not follow the log-logistic distribution, or any distribution for which I can test. Is there a procedure to standardize values that do not follow a distribution? What I am trying to do is calculate the drought index values to use in an attempt to correlate the drought index with irrigation water demand. Thanks in advance. ![histogram](http://i1252.photobucket.com/albums/hh575/jedmoorhead/HistogramofPD1ssCorn.jpg)


title: Treating multiple Dichotomies combined?

text: Let's say I am interested in choosing a new country $c_1, \ldots, c_k$ to live in. For some reason I can only apply to one country and only once.   I know for each country a set of 2000 observations (persons), each consisting of 20 attributes and for each I know whether they were accepted or rejected. (I also know my values for each of those attributes.) My two questions are:   1) I am not sure wether to estimate the probability distribution for acceptance/rejection independently for each country, or not. The motivation for combining observations (or parts of it) of different countries might be that countries probably have similar guidelines for rejection (for example for a convict).   2) Given the training size (2000) and the number of dimensions (20), would a discriminative model (logistic regression,... ) be appropriate or would the amount of data be enough for a generative model like (Naive Bayes). Or would any other approach might be right?    I of course just made this example up to make my question less abstract.


title: Crash course in robust mean estimation

text: I have a bunch (around 1000) of estimates and they are all supposed to be estimates of long-run elasticity. A little more than half of these is estimated using method A and the rest using a method B. Somewhere I read something like "I think method B estimates something _very_ different than method A, because the estimates are much (50-60%) higher". My knowledge of robust statistics is next to nothing, so I only calculated the sample means and medians of both samples... and I immediately saw the difference. Method A is very concentrated, the difference between median and mean is very little, but method B sample varied wildly. I concluded that the outliers and measurement errors skew the method B sample, so I threw away about 50 values (about 15%) that were very inconsistent with theory... and suddenly the means of both samples (including their CI) were very similar. The density plots as well. (In the quest of eliminating outliers, I looked at the range of sample A and removed all sample points in B that fell outside it.) I would like you to tell me where I could find out some basics of robust estimation of means that would allow me to judge this situation more rigorously. And to have some references. I do not need very deep understanding of various techniques, rather read through a comprehensive survey of the methodology of robust estimation. I t-tested for significance of mean difference after removing the outliers and the p-value is 0.0559 (t around 1.9), for the full samples the t stat was around 4.5. But that is not really the point, the means can be a bit different, but they should not differ by 50-60% as stated above. And I don't think they do.


title: Help with multilinear regression model selection

text: The model is from sports. I'm trying to predict the win percentage of each team in the next season based on information available before that season. I have about 240 observations of the numbers of wins in various seasons for various teams, along with the corresponding values of possible predictors. I want to build a multilinear regression model using this data and I want to demonstrate the success or lack of success of whatever model I come up with by showing how my model's predictions compare to those of an "expert" in one or two most recent seasons. In my model, I have 4 potential predictors X1, X2, X3, X4. Two of them, say X1 and X2, are _a priori_ reasonable approximations of the target variable, though obtained with different assumptions, so I expected and found that they are highly correlated. I've read some on multicollinearity, so I understand this causes a variety of problems when one tries to include them in the same model. It would make sense to me that the model which includes X1, X3, and X4 and and the model which includes X2, X3, and X4 would be best, assuming that X3 and X4 are significant, etc. but this also seems a bit unsystematic/arbitrary. I wonder if you could offer some advice on model selection in this situation, possible methods to use, etc. Particularly, I'm confused about what the appropriate way to partition the data is. Obviously if I'm going to be demonstrating the model on a recent season, I should exclude that recent season from the data I train the model on. But should I also set aside some of the data for validation? As you can probably tell, I'm new to modelling. Edit: More details, in light of Rolando's answer: X1 is that team's win percentage the previous season, while X2 is an estimation of next season's winning percentage based on solving the formula for MP/MW for winning percentage and assuming that MP/MW next season will be similar to what it was the previous season. So one estimate is based on an assumption of "continuity" of winning percentage from one season to the next and the other is based an assumption of the "continuity" of this MP/MW quantity from one season to the next, where MP/MW is supposed to be an estimate of the efficiency of the team's organization, basically. Another question: Supposing you have a target variable Y and some number of predictor variables X1, ..., Xn which are already predictions of Y, it's reasonable to think that we can aggregate these predictions to come up with a better estimate of Y. What are some simple methods of doing this which I should read about? Thanks.


title: How to choose asymmetric biplots in correspondence analysis

text: I have categorical data on Police Stations against Crime categories on a 2-way contingency table. I have an asymmetric row map and column map from correspondence analysis using XLSTAT. Which one of the two maps should be used to interpret the results?


title: How to use the Laplace approximation to fit a spatial GLM model

text: I would like to fit a spatial generalized linear mixed model, for example, to model spatial count responses by using the Laplace approximation in available statistical software. I am particularly interested in implementations in R.


title: Autocorrelation testing on non-timeseries data

text: I have a regression model that does not have any time variables and the ordering of the variables and their values does not matter when doing the regression. Are autocorrelation tests like Durbin-Watson still applicable to this kind of data?


title: How to graph interaction effects for panel data

text: I have a pretty generic question which I am guessing could be relevant to many social scientists who deal with panel data sets. What are the best practices for making graphs about interaction effects. There is a variety of articles that talk about packages like ggplot, ggplot2, sjPlot, effects and so on but to my knowledge none of them really takes into account multiple years and multiple individuals. In my field of management what I see most often is the complete disregard for year and individual effects so that an interaction for instance gets fitted by one x-variable on the x-axis, a simple high versus low of the other interacting variable in the graph and the y-variable on the y-axis of course. A difference in slope then shows how the interaction results in different impact on the y-variable. This seems pretty simplistic as it (if I understand correctly) pretty much erases everything idiosyncratic about panel data, so basically it disregards the fundamental reasons why panels are interesting for statistical inference. So hence my question to the broader community: How do you plot results (and interactions) from panel data. I'd appreciate both general and coding suggestions :) Some extra questions:   1. When you have a bunch of control variables, do you chose artificial values for all of them? What about time and individual variables?   2. Is there a way of plotting the difference in residuals between the regression without interaction effects and the ones with interactions (similar to marginal effects) that clarifies how the interacting variables affect the dependent (response) variable y?   3. While uncommon in management, can we draw 3 dimensional figures to show interactions? Looking forward to your feedback and suggestions


title: error variance calculations for reliability analysis of a composite metric in HLM

text: I am trying to determine how to obtain within-group variance for a composite measure based on a set of (weighted) proportions. I have 50 groups being compared on 8 proportion measures, with 4 measures weighted as 1/5 and 4 measures weighted at 1/20. There is a composite score for each group that is the 8 rates multiplied by the weights and summed. Now I need the variance of each group's composite score in order to use in a HLM. The RAND tech report "The Reliability of Provider Profiling" has a nice overview (and SAS code) for individual binomial measures, but the tutorial turns suddenly obscure for assessing the reliability of the composite, particularly with weighted proportions. On p. 31: "For the weighted combination of several binary measures, the sampling variance can be calculated by multiplying the variances of individual measures by the squared weights and summing the components." MAIN QUESTION: So, does that mean a) the individual measure's variance is (p(1-p))/n) since it is a proportion, and b) I'd multiply that by the weight squared, so c) the weighted variances for each measure are summed to make the final errorvariance value for each group's composite score? It seems like that's what the tutorial implies. Or do you need to divide the summed, weighted variance by N to obtain the final variance value? (I've seen that formula in passing.) Or are neither appropriate? Related questions how are variances obtained for these two situations: 1) what if your composite is an O/E where E is the overall mean? 2) what if your data are continuous?


title: Confirmatory factor analysis using SEM: What do we do with items with zero variance?

text: When I run my analysis with items with zero variance, `sem` in `R` would throw up a warning saying that items with zero variance exist. I then move on to remove these from the analysis. My question is:   1. Are there any other ways of resolving this problem besides removing the items perhaps like adding a very small amount of jitter?   2. What literature is there to justify such a procedure? (I have searched but have not come up with anything. I am guessing either this is not a common problem or the solution is too obvious.)


title: Normalizing reading times using analysis of variance

text: I'm trying to compare gaze durations measured using an eye tracker for each word. Since longer words will naturally lead to longer gaze durations, simply comparing the values is inaccurate.               gaze duration    word length     100ms            6     250ms            8     150ms            7      A simple way to normalize the values is by dividing the duration by the number of characters. However, this transformation assumes that reading time is normally a linear increasing function of the number of character, with a value of zero when the number of characters is zero which is not quite right.               gaze duration    word length    Adjusted gaze duration (duration / length)     100ms            6              16.67 ms/character     250ms            8              31.25 ms/character     150ms            7              21.4 ms/character      The best way to normalize the gaze durations was described in this paper: > A more nearly adequate position assumes that reading time is normally a > linear function of number of characters, with a zero intercept. A linear > regression analysis could be used to estimate the slope and zero intercept > of such a function, and thus to estimate the expected reading times for > regions of varying lengths. Deviations from these expected times would > indicate the existence of factors that speeded or slowed the reading of any > given segment. > > Such an analysis was performed by computing the linear regression equation > expressing reading time for each segment in each experimental passage as a > function of the number of characters in it for each subject. The correlation > averaged over all subjects was .38. The regression equation was used to > obtain expected reading time on the basis of number of characters alone for > each segment. The expected reading times were then subtracted from the > obtained reading times and the resulting difference scores were submitted to > an analysis of variance.               linear fit {6,100},{8,250},{7,150} = 75x - 358.333          gaze duration    word length    expected gaze duration    difference     100ms            6              75(6) - 358.333 = 91.66   8.34     250ms            8              75(8) - 358.333 = 241.67  8.33     150ms            7              75(7) - 358.333 = 166.67  -16.67      **My question** How does submitting the difference scores to an analysis of variance finally allow us to get the transformed values? Is merely using the difference scores as the normalized scores considered an adequate normalization? **Note** The best answer should provide an example of the calculation necessary to get the final corrected values, you may use the sample values provided in the question.


title: Average variation points

text: Having a score board                       Points     Day 1   1000     Day 2   1050     Day 3   1010     Day 4   1010     Day 5   1030      How to calculate the average variation points? i.e. How many points, on average, a person gains or loses per day. This must be something really simple. Yet, I can't figure out.


title: Measure difference in answers on survey depending on item response options

text: For my thesis I sent out a questionnaire that has 10 questions. I set defaults: the first seven question had the default set to the desired answer and the last three questions had the default set to the undesired answer. I want to measure whether people followed the default using SPSS. I am thinking I should use a t-test but this is mainly used to measure a difference between populations, not between questions in one survey. Does anyone have an idea how I can measure if my sample was influenced by the default? An example of something I would like to test is: two questions are similar but have a different default, did they answer them the same or did they follow the default.


title: How to avoid and handle survey non-response?

text: One objective of a survey can be to understand the proportion male vs. female users. To the extent a specific gender correlates to particular use cases of a technology, product design, product/feature prioritization and marketing can be differently managed. A survey of existing users will inherently have error bars based on sample size, frequency of sampling, etc. As I understand it, when asking the question 'are you male or female?' as part of a demographic survey where all fields are required, an additional selection bias is introduced into the survey. Specifically, if a survey is sent to 100 users, of which 50 respond, the proportion of responders that are male may not represent the proportion of 100 customers that were initially surveyed. Assuming there is not a existing sample set of known males and known females whom to survey in order to estimate the responder bias, what are some approaches to prevent or to correct for this error?


title: Can a very bad Coefficient of determination ($R^{2}$) not be indicative of model performance?

text: Thanks in advance for the advice. I am trying to build a generalized linear model that has many predictors. The $R^{2}$ value of the model is quite low (.21), but when I use the model to predict against my validation set I am getting very good results. I was under the impression that a low $R^{2}$ value generally means that the predictive power of a model is low. What could be going on here (I am looking for reasons why a model may make good predictions but have a low coefficient of determination)? My training and validation sets have a similar distribution and I believe my validation and training sets to representative of the whole space.