title: Is partial correlation possible after Spearman's correlation?

text: I have two test groups that conducted an online task measuring response times (avg, avg(congruent), avg(incongruent)). I expected one group to be faster than the other but it turned out exactly the opposite way (significant). So now I'm trying to find one or more factors that are responsible for this unexpected outcome. My test data looks like this: participants conducted a Simon task and I measured their response times for each trial. In SPSS I have 3 variables, one for their average response times over all trials, one with their average response time for the trials that were congruent (place on screen and direction match), one with their average response time for the trials that were incongruent (place on screen and direction do not match). My groups are monolingual and multilingual. Literature suggests that multilinguals would score either higher averages overall, or higher averages for just the incongruent trials. My outcomes show that the monolinguals have a higher overall average and that the difference is the same between congruent and incongruent trials for both groups. For each group I know the following variables:   * age   * education level   * accuracy%   * sex   * lurking variables (yes|no)   * average response times (total, congruent trials, incongruent trials)   * first language Each variable in itself seems to be significant, but none of the variables seem to have a significant influence in combination with the original groups. I want to ensure that my test was not faulty and verify the validity of the test method. Maybe one other factor I didn't account for in my hypothesis is causing these results. It would either help me substantiate that my test is correct and that the results are valid, or provide a new angle for future research. I'm not sure what test to use. For my first comparisons I used Spearman's correlation. I was told once that I may not use partial correlations when using Spearman's correlation. Can anyone help me how to proceed from here?


title: Using principal component scores as predictors in mixed-model

text: Would there be any problem with using principal component analysis (e.g. for reduction of dimensionality) so that principal components scores could be used as predictors in a mixed-model? For non-mixed models this strategy is frequently applied (principal component regression) but I am not sure if it is applicable in the context of mixed-models? Please see below a dummy example in R:               library(lme4)     USArrests$score <- prcomp(USArrests[,-1], scale = TRUE)$x[,1]     USArrests$group[1:25)] <- "A"     USArrests$group[26:50] <- "B"     m1 <- lmer(Murder~1+score+(1|group), data=USArrests)     summary(m1)


title: Feature selection and cross-validation

text: I have recently been reading a lot on this site (@Aniko, @Dikran Marsupial, @Erik) and elsewhere about the problem of overfitting occuring with cross validation - (Smialowski et al 2010 Bioinformatics, Hastie, Elements of statistical learning). The suggestion is that _any_ supervised feature selection (using correlation with class labels) performed outside of the model performance estimation using cross validation (or other model estimating method such as bootstrapping) may result in overfitting. This seems unintuitive to me - surely if you select a feature set and then evaluate your model using _only_ the selected features using cross validation, then you are getting an unbiased estimate of generalized model performance on those features (this assumes the sample under study are representive of the populatation)? With this procedure one cannot of course claim an optimal feature set but can one report the performance of the selected feature set on unseen data as valid? I accept that selecting features based on the entire data set may resuts in some data leakage between test and train sets. But if the feature set is static after initial selection, and no other tuning is being done, surely it is valid to report the cross-validated performance metrics? In my case I have 56 features and 259 cases and so #cases > #features. The features are derived from sensor data. Apologies if my question seems derivative but this seems an important point to clarify. **Edit:** On implementing feature selection within cross validation on the data set detailed above (thanks to the answers below), I can confirm that selecting features prior to cross-validation in this data set introduced a _significant_ bias. This bias/overfitting was greatest when doing so for a 3-class formulation, compared to as 2-class formulation. I think the fact that I used stepwise regression for feature selection increased this overfitting; for comparison purposes, on a different but related data set I compared a sequential forward feature selection routine performed prior to cross- validation against results I had previously obtained with feature selection within CV. The results between both methods did not differ dramatically. This may mean that stepwise regression is more prone to overfitting than sequential FS or may be a quirk of this data set.


title: Calculating distances between centroids in an Primer 6 MDS plot

text: I am just wondering if there is a way to calculate the distances between the centroids generated by Primer 6 in a MDS plot. I am trying to find if the centroid distances change beween samples but I can't find a way to obtain values from the program.


title: What is the intuition behind (M)ANCOVA and when/why should one use it?

text: As per my understanding here's the what/why/when of the following hypotheses tests in a crude sense:   * **t-test** : Used when comparing means between two samples   * **ANOVA (one way)** : Used when you have one dependent variable and one independent (i.e., categorical) variable and you wish to analyze the 'means' (i.e., effects) across multiple groups. Simply stated, multi-way t-tests in essence.   * **ANOVA (two way)** : Similar to one-way except you have two independent (i.e., categorical) variables   * **MANOVA** : ANOVA with multiple dependent variables   * **ANCOVA** : ??   * **MANCOVA** : ?? Intuitively, the concepts/intuition behind (M)ANOVA makes sense and I understand when/how to apply it and why is it necessary. I've just overly simplified my understanding about them above. However, I lack the similar intuition behind (M)ANCOVA.


title: what test for analysing level of agreement of users perception on species population trends from ordinal categorical data

text: This is my first post ever so I beg you please forgive me any mistake (I'm also super newbie in stats). I have this dataset that comes from interviews filled by users (divers) about their perception on resources (six shark species) population trends across four decades. Every user was asked to state their perception of population trend for every shark based on five categories: major decrease, decrease, estable, increase, and major increase. They were asked to constrain their answers per decade: trend for species X within decade 80s, trend for the same species within 90s, 00s, etc for the rest species. Worth to note that the amount of answers per decade varied according to the level of experience of divers, which means that there were more divers with experience toward the last two decades, and just few since the 80s (a total of four divers with experience since 80s; 14 since 90s; and, 25 since 2000s). I have attached a figure on the answers of those 25 divers on one shark species. Obtained scores from interviewing divers about an X shark species population trend across the last four decades. Answers obtained from 25 divers. ![enter image description here](http://i.stack.imgur.com/z1Jfu.png) I've been looking for the simplest way to assess their level of agreement in any specific trend across decades. Find out that ANOVAs are no near a good way to assess any type of categorical data and literature redirected me into logic mixed models, although I still don't have clear the picture how those will let me assess level of agreement and overall species trends across time. Could anybody please give some light to my analysis and tell me which should be the best statistical test to assess level of agreement and trends across time of this dataset? Will logic mixed models work for this? Please, I use R so I would like to kindly request your help and point me in any specific R routine/package. Also, if there is anything I should add to make my question clearer, please do let me know. Thanks so much!! C.


title: Translating R lme comand to mathematical equation

text: I would appreciate if someone could help me in translating the following R command into a mathematical equation:               lme(score ~ factor(timeslot), random=~1|subjectid, data=a)      For each subject there are observations at 6 different times (timeslots) during a day.


title: Prediction method where predictors and response variable are binary

text: I have a group of binary tasks performed by multiple subjects. Every task can be either performed right or wrong (i.e.,1/0). My goal is to predict the accuracy of future task given the performance on previous tasks. As tasks considered to be independent, I thought that it makes sense to average the accuracy by making predictions for every task given the rest. Could you please advise me what is a suitable regression or any predictive model for this kind of data Additionally, I am interested in cross validation between variables: Lets say I have 10 variables. I'm interested in averaged accuracy of prediction when first 9 variables are predictors and the last one is the response variable, then variables 1-8,10 are predictors and 9th is the response variable, and so on such that every variable calculated once as a response variable. Thank you for your help


title: How to select the final model with elastic net feature selection, cross validation and SVM?

text: I have a dataset of some 100 samples, each with >10,000 features, some of which highly correlated. Here's what I am doing currently.   1. Split the data set into three folds.   2. For each fold,   2.1 Run elastic net for 100 values of lambda. (this returns a nfeatures x 100 matrix)   2.2 Take a union of all non-zero weights. (returning a nfeatures x 1 vector)   3. Select features corresponding to the non-zero weights returned from 2.2   4. Use these features for training and testing SVM. My problem is that in step 3, for each fold I get a different set of features. How do I get one final model out of this? One final list of relevant features? Can I take an intersection of the selected features in step 3 for all folds? Features that are selected in all three folds would appear to be the most stable/significant. Can I do this, or is it cheating?


title: Appropriate accuracy test if drawing distributions on a set of actual values and predicted values?

text: My question regards statistics as applied to text mining. I have used substring matching to determine the predicted set of keywords. I then classify those keywords into broader groups. Given a body of text, I want to draw a single group based on the keywords. I used the Multinomial-Dirichlet to do so, where the Dirichlet is updated based on the counts of a certain type of keyword. First I draw the dirichlet probability vector and then use this to make the multinomial draw. Repeating this pair of drawing several thousand times, I examine the average proportion found in each category. I do the same to find average proportion for the actual known set of keywords. To compare my results, which statistical test would be best suited for accuracy/error analysis?


title: Multiple values in a timestamp

text: I am working with a time series (discrete) that has ideally 1 value per time stamp. In some cases, we are seeing multiples that have a wide range all recorded with the same time stamp. Up to now, we had been averaging the data but have found quite a few instances of very big ranges. Should we instead be considering the modes (sometimes there are only 2 values), should we be considering a different method? When should we consider just throwing them out? There are many missing values in this data so using previous and following values in the sequence is not ideal. I can't recall from my previous stats what the preferred method of dealing with these values is, as I haven't come across this in some time (maybe ever, in practice). Any ideas on how to handle the situation is greatly appreciated. Resources to look up would very much help as well! Essentially the data sometimes looks like this, where the values are separated by ";" (these are not in sequence):                00:43:00      "78;66;81"      02:01:00      "68;74"      03:53:00      "78;86;88;95;95;111;102;97;94;95"      15:48:00      "81;120;97;58;84      Edit: While the ideal is to have 1 per timestamp, the reason why there is more than one in some instances is the truncation at minutes by the database. However, we need to pick just one value to represent that timestamp as we don't care about the seconds details.


title: AUC in ordinal logistic regression

text: I'm using 2 kind of logistic regression - one is the simple type, for binary classification, and the other is ordinal logistic regression. For calculating the accuracy of the first, I used cross-validation, where I computed the AUC for each fold and than calculated the mean AUC. How can I do it for the ordinal logistic regression? I've heard about generalized ROC for multi-class predictors, but I'm not sure how to compute it. Thanks!


title: How to correlate two time series, with possible time differences

text: The question is similar to How to correlate two time series with gaps and different time bases? but with regular sample frequency and identical time base between two time series. Say I have two time series, one is the stock price of a given company, the other is the ratio of positive comment about the company on social networks. How can I statistically verify they are correlated? Assume they are both sampled weekly or daily. The challenge lies in the fact there might be time difference from each other. In the context of my example above, the positive evaluation from consumer might be later reflect on the stock price. There are two things I want to investigate in this problem: 1) is the two time series are correlated? 2) what is the time window to maximize their correlation?


title: Correlation assumption of Hochberg procedure in repeated measures design

text: I'm running a 3-factor (Between by Between by Within) ANOVA with the correlation structure modeled using a Random Intercept and Gaussian Serial Correlation. I have the following variables: Dependent variable: Bodyweight Independent: Size - size of the subject (between). 2 levels: normal and large   Diet - diet given (between). 2 levels: normal or experimental diet.   Time - weekly bodyweight measurements (within). 12 levels: 1:12 I have a significant 3-way interaction and all 2-way interactions(p<0.05), and so I am stratifying the 3-factors to test all pairwise differences between Normal and Large Size, on each Diet, and at each week. I'm doing these specific comparisons because I'm interested in whether Size further increases the weight gain known to occur on the experimental diet. Using this approach, I have a total of 24 post-hoc comparisons. I would like to improve power by using Hochberg's procedure rather than Holm's to adjust for multiple comparisons, but I'm unclear whether my data meet the assumption of independent or positively associated test statistics/p-values (as per Simes' test and multcomp documentation in R). I'm employing the mixed model because of the positively correlated repeated measures, so I'm tempted to conclude that the post-hoc/pairwise test statistics at each time would also be positively correlated. But I could very well be wrong. Here are my questions:   1. Would the pairwise test statistics/p-values of these repeated measures data meet the independence/positive assumption of Hochberg's procedure?   2. Is there a way to test or at least calculate/extract the correlations of the test statistics/p-values to empirically validate the assumption?


title: Correlation among categories of categorical nominal variables

text: I have a data set with two categorical nominal variables (both with 5 categories). I would like to know if (and how) I am able to identify potential correlations among categories within these variables. In other words whether for example the results of category 1 in variable 1 show a strong correlation with a specific category in variable 2. Since I have two variables with 5 categories, the total correlation analysis for all the categories would come down to 25 results (at least if it works the way I hope/expect it to work). I have tried to formulate the problem into concrete questions: Question 1: Let's say I transfer the categorical variable into 5 different dummy variables per value. This same procedure I run for the second variable as well. Then I want to determine the correlation between dummy 1.1 and 2.1 (for example). Is it statistically correct for me to execute this procedure by means of an ordinary correlation coefficient procedure? Does the correlation coefficient resulting from this procedure provide a proper insight in a correlation among the two dummy variables? Question 2: If the procedure described in question one is a valid procedure, is there a way to execute this analysis for all categories of 2 (or perhaps more) categorical nominal variables all at once? The program I am using is SPSS (20).


title: Which statistical test should I use to test for associations between phenotype and genotype?

text: I need to select a statistical test in order to answer the question of whether the `TASR238` genotype is a good predictor of the ability to taste PTC (phenotype) Here are the data I collected for my lab report:                   PTC taster           Genoty   Strong  Weak   Non     T/T       10     4      2     T/t       30     26     12     t/t       5      4      18      (also found here) I think I have to use $\chi^2$, but would it be an $\chi^2$ test for differences? How do I calculate expected phenotype based on genotype, would it just be the total number of genotype?


title: How to efficiently expand a simple random sample after population growth?

text: I have a simple random sample (SRS) $S$ of size $n$ drawn from a population $D_1$ containing $N_1$ members. A new set $D_2$ with $N_2$ members is added to the population. I would like to create a new set $S'$ that is a SRS from $D_1 \cup D_2$. One strategy that won't work is to draw a SRS of size $nN_2/N_1$ (or any other size) from $D_2$ and add is to $S$. (The result would be a stratified random sample from $D1 \cup D2$, but not a SRS from it.) A strategy that will work is to first make one trial on a random variable with a hypergeometric distribution (population size $N_1+N_2$, number of draws $n'$, and number of success states in the population $N_1$). The resulting value $n'_{1}$ is the number of items in $S'$ that should be taken from $D_1$. If $n'_{1} < n$, we discard a SRS of size $n - n'_{1}$ from $S$ to form the initial $S'$. If $n'_{1} > n$, we select an additional SRS of size $n'_{1} - n$ from $D_1 \setminus S$, and put these items, plus all items from $S$ in $S'$. We also draw an SRS of size $n'_{2} = n' \- n'_{1}$ from $D_2$ and add it to $S'$. I'm interested in choosing $n'$ to minimize $E[|S' \setminus S|]$, the expected size of $S' \setminus S$. In other words, the expected value of the number of items from $D_2$ plus the number of new items, if any, from $D_1$. I have three questions: Q1. If we have the additional constraint that $n' \ge j$ for some $j$, is there a simple argument that the optimal $n'$ is $j$? Seems like it has to be, but I'd rather not work out the expected value formula and slog through the derivatives. Q2. Does anyone know of a better algorithm for the above problem? Q3. The problem I'd _really_ like to solve is more complex. Each item in $D_1$ and $D_2$ is either an A or a B. $S$ has at least $k$ A's, and I want to minimize $E[|S' \setminus S|]$ subject to the constraint that $S'$ has at least $k$ A's. The above algorithm no longer works, because we can't pick $n'$ in advance. (We don't know how rich $D_2$ is in A's, and we might get unlucky anyway.) Is there a simple algorithm for this case? I've perused Thompson, S. K. and Seber, G. A. G. _Adaptive Sampling_. Wiley, New York, 1996 but can't find exactly this problem treated. But my gut tells me that both the above problems have been faced and solved many times before.


title: Binomial regression asymptotes

text: Binomial logistic regression has upper and lower asymptotes of 1 and 0 respectively. However, accuracy data (just as an example) may have upper and lower asymptotes vastly different to 1 and/or 0. I can see three potential solutions to this:   1. Don't worry about it if you are getting good fits within the area of interest. If you aren't getting good fits then:   2. Transform the data so that the minimum and maximum number of correct responses in the sample give proportions of 0 and 1 (instead of say 0 and 0.15).   _or_   3. Use non-linear regression so that you can either specify the asymptotes or have the fitter do it for you. It seems to me that options 1 & 2 would be preferred over option 3 largely for simplicity reasons, in which case option 3 is perhaps the better option because it can yield more information? _edit_   Here's an example. Total possible correct for accuracy is 100, but the maximum accuracy in this case is ~ 15.               accuracy <- c(0,0,0,0,0,1,3,5,9,13,14,15,14,15,16,15,14,14,15)     x<-1:length(accuracy)     glmx<-glm(cbind(accuracy, 100-accuracy) ~ x, family=binomial)     ndf<- data.frame(x=x)     ndf$fit<-predict(glmx, newdata=ndf, type="response")     plot(accuracy/100 ~ x)     with(ndf, lines(fit ~ x))      Option 2 (as per comments and to clarify my meaning) would then be the model               glmx2<-glm(cbind(accuracy, 16-accuracy) ~ x, family=binomial)      Option 3 (for completeness) would be something akin to:               fitnls<-nls(accuracy ~ upAsym + (y0 - upAsym)/(1 + (x/midPoint)^slope),        start = list("upAsym" = max(accuracy), "y0" = 0, "midPoint" = 10, "slope" = 5),        lower = list("upAsym" = 0, "y0" = 0, "midPoint" = 1, "slope" = 0),        upper = list("upAsym" = 100, "y0" = 0, "midPoint" = 19, hillslope = Inf),        control = nls.control(warnOnly = TRUE, maxiter=1000),       algorithm = "port")


title: Multiple linear regression question

text: I am running a multiple regression of the form _Y_ = $\beta_0$ + $\beta_1$*$X_1$ + $\beta_2$*$X^2_1$ + $\beta_3$*$X_2$ + $\beta_4$*$X_3$ on a time-series dataset. I want to plot the relationship between _Y_ and $X_1$ based on this equation such that I get a curve after taking into account the effects of $X_2$ and $X_3$ (or in other words, after the effects of $X_2$ and $X_3$ are removed from _Y_ ). What is the most appropriate way to do this?


title: 5-sigma anomaly problem

text: Struggling with how to start this question: A factory processes chickens with masses between 1.5 and 2.5kg. "Natural" chickens have masses equally spread across this range. However, "GM" chickens have masses between 2.3 and 2.4kg. If a factory produces 10 million chickens per year, how many GM chickens would be required to produce a 5-sigma anomaly in the chicken mass spectrum? I understand that the chicken masses are going to follow a normal distribution but quite lost after that. Anybody want to shed some light?


title: Two-dimensional Kolmogorov-Smirnov

text: I would like to run some two dimensional Kolmogorov-Smironov tests to determine whether a two-dimensional distribution fits with a reference. Is there any package or application that I could I use in a relatively straightforward fashion? Or is there a different algorithm that is preferrable? I have just a basic statistical knowledge.


title: Test for martingale

text: Given a time series, I am trying to determine whether the underlying stochastic process is a martingale or not. Until now, I have deduced that since I know that the process is always bounded after a deterministic time, if I can prove that it is a Brownian motion, I will be done. Is there a standard test for this?


title: What is the intuition for dependence assumption in Benjamini and Hochberg (1995)?

text: For dealing with multiple testing, Benjamini and Hochberg (1995) assume that the tests are independent. Benjamini and Yekutieli (2001) show that (1) the 1995 result also holds under (I think) weaker conditions and that (2) including the factor $\Sigma^m_{i=1}\frac{1}{i}$ makes the procedure applicable for any dependence structure (see Theorem 1.3). Since this summation is greater than or equal to 1, it would make the procedure more conservative. My question is why is independent tests not the worst case assumption? Is there some simple example to illustrate why not? Does it involve a negative dependence structure of some sort? Maybe my understanding is faulty, but I thought that the Bonferroni adjustment (for the familywise error) worked for any dependence and that independence is the worst case assumption. Perhaps this is not correct and hence my confusion with BH.


title: How exactly can the Cox-model ignore exact times?

text: The Cox model does not depend on the times itself, instead it only needs an ordering of the events. How come it doesn't need the time, as all of the models I've seen so far are dependent on the exact time point and/or interval?


title: Fitting mixture distributions and computing goodness-of-fit?

text: This question is a follow-up from a previous question of mine here. Thanks to @Glen_b, @gung and @rbatt for teaching me so many new things yesterday. It was mentioned in passing that mixture distributions might be a possibility so I figured I will learn about these. I have seen some articles on how to fit mixture distributions to data and they seem very interesting. To fit two lognormals to my data, I did the following:               library(mixtools)          x <- c(1528L, 285L, 87138L, 302L, 115L, 416L, 8940L, 19438L, 165820L, 540L, 1653L, 1527L, 974L, 12999L, 226L, 190L, 306L, 189L, 138542L, 3049L, 129067L, 21806L, 456L, 22745L, 198L, 44568L, 29355L, 17163L, 294L, 4218L, 3672L, 10100L, 290L, 8341L, 128L, 11263L, 1495243L, 1699L, 247L, 249L, 300L, 351L, 608L, 186684L, 524026L, 1392L, 396L, 298L, 1063L, 11102L, 6684L, 6546L, 289L, 465L, 261L, 175L, 356L, 61652L, 236L, 74795L, 64982L, 294L, 95221L, 322L, 38892L, 2146L, 59347L, 2118L, 310801L, 277964L, 205679L, 5980L, 66102L, 36495L, 580277L, 27600L, 509L, 21795L, 21795L, 301L, 617L, 331L, 250L, 123501L, 144L, 347L, 121443L, 211L, 232L, 445783L, 9715L, 10308L, 1921L, 178L, 168L, 291L, 6915L, 6735L, 1008478L, 274L, 20L, 3287L, 591208L, 797L, 586L, 170613L, 938L, 3121L, 249L, 1497L, 24L, 1407L, 1217L, 1323L, 272L, 443L, 49466L, 323L, 323L, 784L, 900L, 26814L, 2452L, 214713L, 3668L, 325L, 20439L, 12304L, 261L, 137L, 379L, 2273L, 274L, 17760L, 920699L, 13L, 485644L, 1243L, 226L, 20388L, 584L, 17695L, 1477L, 242L, 280L, 253L, 17964L, 7073L, 308L, 260692L, 155L, 58136L, 16644L, 29353L, 543L, 276L, 2328L, 254L, 1392L, 272L, 480L, 219L, 60L, 2285L, 2676L, 256L, 234L, 1240L, 219714L, 102174L, 258L, 266L, 33043L, 530L, 6334L, 94047L, 293L, 536L, 48557L, 4141L, 39079L, 23259L, 2235L, 17673L, 28268L, 112L, 64824L, 127992L, 5291L, 51693L, 762L, 1070735L, 179L, 189L, 157L, 157L, 122L, 1045L, 1317L, 186L, 57901L, 456126L, 674L, 2375L, 1782L, 257L, 23L, 248L, 216L, 114L, 11662L, 107890L, 203022L, 513L, 2549L, 146L, 53331L, 1690L, 10752L, 1648611L, 148L, 611L, 198L, 443L, 10061L, 720L, 10L, 24L, 220L, 38L, 453L, 10066L, 115774L, 97713L, 7234L, 773L, 90154L, 151L, 1560L, 222L, 51558L, 214L, 948L, 208L, 1127L, 221L, 169L, 1528L, 78959L, 61566L, 88049L, 780L, 6196L, 633L, 214L, 2547L, 19088L, 119L, 561L, 112L, 17557L, 101086L, 244L, 257L, 94483L, 6189L, 236L, 248L, 966L, 117L, 333L, 278L, 553L, 568L, 356L, 731L, 25258L, 127931L, 7735L, 112717L, 395L, 12960L, 11383L, 16L, 229067L, 259076L, 311L, 366L, 2696L, 7265L, 259076L, 3551L, 7782L, 4256L, 87121L, 4971L, 4706L, 245L, 34457L, 4971L, 4706L, 245L, 34457L, 258L, 36071L, 301L, 2214L, 2231L, 247L, 537L, 301L, 2214L, 230L, 1076L, 1881L, 266L, 4371L, 88304L, 50056L, 50056L, 232L, 186336L, 48200L, 112L, 48200L, 48200L, 6236L, 82158L, 6236L, 82158L, 1331L, 713L, 89106L, 46315L, 220L, 5634L, 170601L, 588L, 1063L, 2282L, 247L, 804L, 125L, 5507L, 1271L, 2567L, 441L, 6623L, 64781L, 1545L, 240L, 2921L, 777L, 697L, 2018L, 24064L, 199L, 183L, 297L, 9010L, 16304L, 930L, 6522L, 5717L, 17L, 20L, 364418L, 58246L, 7976L, 304L, 4814L, 307L, 487L, 292016L, 6972L, 15L, 40922L, 471L, 2342L, 2248L, 23L, 2434L, 23342L, 807L, 21L, 345568L, 324L, 188L, 184L, 191L, 188L, 198L, 195L, 187L, 185L, 33968L, 1375L, 121L, 56872L, 35970L, 929L, 151L, 5526L, 156L, 2687L, 4870L, 26939L, 180L, 14623L, 265L, 261L, 30501L, 5435L, 9849L, 5496L, 1753L, 847L, 265L, 280L, 1840L, 1107L, 2174L, 18907L, 14762L, 3450L, 9648L, 1080L, 45L, 6453L, 136351L, 521L, 715L, 668L, 14550L, 1381L, 13294L, 13100L, 6354L, 6319L, 84837L, 84726L, 84702L, 2126L, 36L, 572L, 1448L, 215L, 12L, 7105L, 758L, 4694L, 29369L, 7579L, 709L, 121L, 781L, 1391L, 2166L, 160403L, 674L, 1933L, 320L, 1628L, 2346L, 2955L, 204852L, 206277L, 2408L, 2162L, 312L, 280L, 243L, 84050L, 830L, 290L, 10490L, 119392L, 182960L, 261791L, 92L, 415L, 144L, 2006L, 1172L, 1886L, 233L, 36123L, 7855L, 554L, 234L, 2292L, 21L, 132L, 142L, 3848L, 3847L, 3965L, 3431L, 2465L, 1717L, 3952L, 854L, 854L, 834L, 14608L, 172L, 7885L, 75303L, 535L, 443347L, 5478L, 782L, 9066L, 6733L, 568L, 611L, 533L, 1022L, 334L, 21628L, 295362L, 34L, 486L, 279L, 2530L, 504L, 525L, 367L, 293L, 258L, 1854L, 209L, 152L, 1139L, 398L, 3275L, 284178L, 284127L, 826L, 751L, 1814L, 398L, 1517L, 255L, 13745L, 43L, 1463L, 385L, 64L, 5279L, 885L, 1193L, 190L, 451L, 1093L, 322L, 453L, 680L, 452L, 677L, 295L, 120L, 12184L, 250L, 1165L, 476L, 211L, 4437L, 7310L, 778L, 260L, 855L, 353L, 97L, 34L, 87L, 137L, 101L, 416L, 130L, 148L, 832L, 187L, 291L, 4050L, 14569L, 271L, 1968L, 6553L, 2535L, 227L, 202L, 647L, 266L, 2681L, 106L, 158L, 257L, 234L, 1726L, 34L, 465L, 436L, 245L, 245L, 2790L, 104L, 1283L, 44416L, 142L, 13617L, 232L, 171L, 221L, 719L, 176L, 5838L, 37488L, 12214L, 3780L, 5556L, 5368L, 106L, 246L, 101L, 158L, 10743L, 5L, 46478L, 5286L, 9866L, 32593L, 174L, 298L, 19617L, 19350L, 230L, 78449L, 78414L, 78413L, 78413L, 6260L, 6260L, 209L, 2552L, 522L, 178L, 140L, 173046L, 299L, 265L, 132360L, 132252L, 4821L, 4755L, 197L, 567L, 113L, 30314L, 7006L, 10L, 30L, 55281L, 8263L, 8244L, 8142L, 568L, 1592L, 1750L, 628L, 60304L, 212553L, 51393L, 222L, 13471L, 3423L, 306L, 325L, 2650L, 74796L, 37807L, 103751L, 6924L, 6727L, 667L, 657L, 752L, 546L, 1860L, 230L, 217L, 1422L, 347L, 341055L, 4510L, 4398L, 179670L, 796L, 1210L, 2579L, 250L, 273L, 407L, 192049L, 236L, 96084L, 5808L, 7546L, 10646L, 197L, 188L, 19L, 167877L, 200509L, 429L, 632L, 495L, 471L, 2578L, 251L, 198L, 175L, 19161L, 289L, 20718L, 201L, 937L, 283L, 4829L, 4776L, 5949L, 856907L, 2747L, 2761L, 3150L, 3142L, 68031L, 187666L, 255211L, 255231L, 6581L, 392991L, 858L, 115L, 141L, 85629L, 125433L, 6850L, 6684L, 23L, 529L, 562L, 216L, 1450L, 838L, 3335L, 1446L, 178L, 130101L, 239L, 1838L, 286L, 289L, 68974L, 757L, 764L, 218L, 207L, 3485L, 16597L, 236L, 1387L, 2121L, 2122L, 957L, 199899L, 409803L, 367877L, 1650L, 116710L, 5662L, 12497L, 613889L, 10182L, 260L, 9654L, 422947L, 294L, 284L, 996L, 1444L, 2373L, 308L, 1522L, 288L, 937L, 291L, 93L, 17629L, 5151L, 184L, 161L, 3273L, 1090L, 179840L, 1294L, 922L, 826L, 725L, 252L, 715L, 6116L, 259L, 6171L, 198L, 5610L, 5679L, 862L, 332L, 1324L, 536L, 98737L, 316L, 5608L, 5526L, 404L, 255L, 251L, 14067L, 3360L, 3623L, 8920L, 288L, 447L, 453L, 1604687L, 115L, 127L, 127L, 2398L, 2396L, 2396L, 2398L, 2396L, 2397L, 154L, 154L, 154L, 154L, 887L, 636L, 227L, 227L, 354L, 7150L, 30227L, 546013L, 545979L, 251L, 171647L, 252L, 583L, 593L, 10222L, 2660L, 1864L, 2884L, 1577L, 1304L, 337L, 2642L, 2462L, 280L, 284L, 3463L, 288L, 288L, 540L, 287L, 526L, 721L, 1015L, 74071L, 6338L, 1590L, 582L, 765L, 291L, 983L, 158L, 625L, 581L, 350L, 6896L, 13567L, 20261L, 4781L, 1025L, 722L, 721L, 1618L, 1799L, 987L, 6373L, 733L, 5648L, 987L, 1010L, 985L, 920L, 920L, 4696L, 1154L, 1132L, 927L, 4546L, 692L, 702L, 301L, 305L, 316L, 313L, 801L, 788L, 14624L, 14624L, 9778L, 9778L, 9778L, 9778L, 757L, 275L, 1480L, 610L, 68495L, 1152L, 1155L, 323L, 312L, 303L, 298L, 1641L, 1607L, 1645L, 616L, 1002L, 1034L, 1022L, 1030L, 1030L, 1027L, 1027L, 934L, 960L, 47L, 44L, 1935L, 1925L, 43L, 47L, 1933L, 1898L, 938L, 830L, 286L, 287L, 807L, 807L, 741L, 628L, 482L, 500L, 480L, 431L, 287L, 298L, 227L, 968L, 961L, 943L, 932L, 704L, 420L, 548L, 3612L, 1723L, 780L, 337L, 780L, 527L, 528L, 499L, 679L, 308L, 1104L, 314L, 1607L, 990L, 1156L, 562L, 299L, 16L, 20L, 287L, 581L, 1710L, 1859L, 988L, 962L, 834L, 1138L, 363L, 294L, 2678L, 362L, 539L, 295L, 996L, 977L, 988L, 39L, 762L, 579L, 595L, 405L, 1001L, 1002L, 555L, 1102L, 54L, 1283L, 347L, 1384L, 603L, 307L, 306L, 302L, 302L, 288L, 288L, 286L, 292L, 529L, 56844L, 1986L, 503L, 751L, 3977L, 367L, 4817L, 4631L, 4609L, 4579L, 937L, 402L, 257L, 570L, 1156L, 3297L, 3948L, 4527L, 3119L, 15227L, 3893L, 538L, 802L, 5128L, 595L, 522L, 1346L, 449L, 443L, 323L, 372L, 369L, 307L, 246L, 260L, 342L, 283L, 963L, 751L, 108L, 280L, 320L, 287L, 285L, 283L, 529L, 536L, 298L, 29427L, 29413L, 761L, 249L, 255L, 304L, 297L, 256L, 119L, 288L, 564L, 234L, 226L, 530L, 766L, 223L, 5858L, 5568L, 481L, 462L, 8692L, 498L, 330L, 7604L, 15L, 121738L, 121833L, 826L, 760L, 208937L, 1598L, 1166L, 446L, 85598L, 513L, 84897L, 50239L, 308L, 1351L, 283L, 7100L, 7101L, 321L, 1019L, 287L, 253L, 634L, 629L, 628L, 678L, 1391L, 1147L, 853L, 287L, 1174L, 287L, 197145L, 197116L, 147L, 147L, 712L, 274L, 283L, 907L, 434L, 1164L, 30L, 599L, 577L, 315L, 1423L, 1250L, 30L, 1502L, 296L, 348L, 617L, 339L, 328L, 123L, 338L, 332L, 47133L, 288L, 340L, 1524L, 1049L, 1072L, 1031L, 1059L, 1038L, 989L, 52L, 54L, 986L, 46L, 1202L, 1272L, 43L, 785L, 761L, 16924L, 289L, 264L, 453L, 365L, 356L, 280L, 16520L, 281L, 255L, 244L, 642L, 1003L, 951L, 921L, 1011L, 45L, 932L, 973L, 39L, 40L, 159L, 566L, 49L, 1161L, 50L, 200L, 215L, 361L, 377L, 980L, 935L, 882L, 281L, 280L, 1025L, 319L, 690L, 284L, 271L, 276L, 286L, 371L, 324L, 304L, 311L, 341L, 603L, 11566L, 270L, 286L, 342L, 326L, 11018L, 282L, 271L, 286L, 586L, 604L, 750L, 608L, 523L, 506L, 3303L, 1079797L, 1079811L, 530L, 2631L, 882L, 628L, 30L, 11905L, 12966L, 390995L, 322353L, 1763L, 1755L, 709L, 713L, 365L, 351L, 205L, 393L, 284L, 39417L, 320L, 322L, 8039L, 995L, 625L, 785L, 298L, 518L, 467L, 1050L, 329L, 141345L, 55566L, 40318L, 287L, 220L, 309346L, 220L, 215314L, 304L, 296L, 4301L, 4311L, 1543L, 1549L, 2876L, 2894L, 287L, 290L, 215L, 605L, 577L, 254L, 1330L, 1863L, 140L, 328L, 284L, 291L, 283L, 1701L, 1696L, 519L, 499L, 2440007L, 289L, 294L, 311L, 324L, 4793L, 4808L, 249L, 205L, 219L, 638L, 2653L, 2648L, 351L, 323L, 1056L, 327L, 794L, 1491L, 284L, 289L, 220L, 765L, 565L, 808L, 832L, 772L, 41668L, 42307L, 6843L, 6612L, 6598L, 241164L, 531L, 554L, 1246L, 459L, 971504L, 805L, 2615L, 2290L, 2086L, 2063L, 2685L, 2704L, 275L, 461L, 458L, 317L, 889L, 335L, 974L, 959L, 253142L, 257L, 250L, 282L, 293L, 666L, 4991L, 287L, 588L, 555L, 3585L, 3195L, 481L, 2405L, 135266L, 571L, 1805L, 365L, 340L, 232L, 224L, 298L, 3682L, 3677L, 577L, 571L, 288L, 297L, 293L, 291L, 256L, 214L, 1257L, 1271L, 65471L, 65471L, 65476L, 65476L, 4680L, 4675L, 339L, 329L, 284L, 288L, 4859L, 4851L, 2534L, 24222L, 330684L, 330684L, 2116L, 282L, 412L, 429L, 2324L, 1978L, 502L, 286L, 943149L, 256L, 288L, 286L, 1098L, 1125L, 442L, 240L, 182L, 2617L, 1068L, 25204L, 170L, 418L, 1867L, 8989L, 1804L, 1240L, 6610L, 1237L, 1750L, 1565L, 1565L, 3662L, 1803L, 218L, 172L, 780L, 1418L, 2390L, 7514L, 23214L, 1464L, 1060L, 1503L, 308802L, 308357L, 21691L, 298817L, 289875L, 4442L, 289284L, 235L, 456L, 676L, 897L, 289109L, 1865L, 288030L, 287899L, 287767L, 287635L, 286639L, 286509L, 286157L, 1427L, 2958L, 4340L, 5646L, 282469L, 7016L, 279353L, 278568L, 316L, 558L, 3501L, 1630L, 278443L, 1360L, 828L, 1089L, 278430L, 278299L, 278169L, 278035L, 277671L, 277541L, 277400L, 277277L, 276567L, 285L, 555L, 834L, 1084L, 1355L, 5249L, 14776L, 1441L, 755L, 755L, 70418L, 3135L, 1026L, 1497L, 949663L, 68L, 526058L, 1692L, 150L, 48370L, 4207L, 4088L, 197551L, 197109L, 196891L, 196634L, 2960L, 194319L, 194037L, 3008L, 3927L, 178762L, 178567L, 403L, 178124L, 2590L, 177405L, 177179L, 301L, 328L, 390685L, 390683L, 575L, 1049L, 819L, 367L, 289L, 277L, 390L, 301L, 318L, 3806L, 3778L, 3699L, 3691L)          mixmdl = normalmixEM(log(x), k=2, epsilon = 1e-08, maxit = 1000, maxrestarts=20,                           verb = TRUE, fast=FALSE, ECM = FALSE, arbmean = TRUE, arbvar = TRUE)          plot(mixmdl, which=2))     lines(density(log(x)), lty=2,lwd=2)      And, the result looks like this, which seems decent enough visually: ![enter image description here](http://i.stack.imgur.com/Quxfc.png) I was wondering how I can test the goodness-of-fit for this mixture. I was looking into the approach outlined here by jbowman but I am unable to get it to work. I changed `pnorm` to `plnorm` because I am testing `lognormal` here in my case.               pmlnorm <- function(x, meanlog, sdlog, pmix) {       pmix[1]*plnorm(x,meanlog[1],sdlog[1]) + (1-pmix[1])*plnorm(x,meanlog[2],sdlog[2])     }          ks.test(log(x), pmlnorm, meanlog=mixmdl$mu, sdlog=mixmdl$sigma, pmix=mixmdl$lambda)          One-sample Kolmogorov-Smirnov test          data:  log(x)      D = 0.9987, p-value < 2.2e-16     alternative hypothesis: two-sided           Warning message:     In ks.test(log(x), pmlnorm, meanlog = mixmdl$mu, sdlog = mixmdl$sigma,  :       ties should not be present for the Kolmogorov-Smirnov test      I think I am doing something fundamentally wrong or maybe these are the actual values. Can someone help me out?


title: I want to learn about ROC curve -- what is the canonical textbook?

text: I want to learn about Receiver-Operator-Characteristic curves, and metrics. I have read through online webpages with some basics, and I have used MATLAB built-ins to create ROC plots. It tells me about area sums. I have also looked at the JMP module for classification and for Logistic plot which make ROC plots. What is the canon here? What is a great textbook with both deep fundamentals and worked out theory? MATLAB and R are preferred. Who made the field? When? Why? I'm looking for good references. Thanks.


title: Good (2d) visualization of a mixture model clustering

text: I have a specific problem which I'm surprised I don't find answers on-line and I hope somebody here has a good suggestion for me. I'm working with a large data set which I'm clustering into specific groups using custom cluster densities. So the original space is quite heterogeneous and has a large number of feature. During and after model optimization, the responsibility vector of each data point is projection of the original features into k dimensions (k = number of clusters). Each such vector also sums to one, which makes them a (discrete) distribution over the space of clusters. Using these soft assignments, I can always create a heatmap visualization of the corresponding confusion matrix (when I have categorical labels). Typically, these assignments are quite sharp, meaning that many vectors are of the form (0,0,1,0,0...) Now, I would like to see a (2d) representation of the points grouped by the information in the responsibility matrix. This way, I could visualize each step in the optimization procedure without knowing actual labels. My first thoughts were PCA, MDS and graph layout algorithms. However, although the heatmap suggests a clear grouping, the PCA does look quite dense with points aligned in lines. As MDS and graphs use distances, I thought about calculating pairwise Hellinger or earth mover distances between the responsibility vectors to apply any of the MDS or graph layout algorithms. However, MDS wasn't successful so far. Has anybody done something like this before? Ideally, I would like to see an animation of the clustering procedure as data points being grouped together. Thanks for your input.


title: Doing a simple Growth Mixture- and Latent Growth Class analysis in R

text: I'm attempting to run a LCGA and GMM analysis in R, but I'm not sure how to go about it. The data I have is pretty simple, so let me explain that first. The dataset is provided by a number of patients who have given item responses (on a Likert scale from 0 to 4) to some questionnaire item at five different time points. So our data is given by               Time  ID    Response     --------------------     1    1       2     2    1       1     3    1       2     4    1       1     5    1       3          1    2       4     2    2       1     3    2       3     4    2       2     5    2       1      and so on, for a total of $n = 64$ individials. I've discovered the `lcmm` package, but I am not sure how to properly use it. Specifically, I am wondering how to properly set the arguments of the `lcmm()` function. Am I correct if I say that I'm doing a Latent Class Growth analysis if I run `lcgm(fixed = Response ~ Time, ng = k, subject = "ID")` and a Growth Mixture analysis if I run `lcgm(random = Response ~ Time, ng = k, subject = "ID")` for a given number of classes $k$? I'm a bit confused because the `fixed` argument seems to be mandatory, but then I'm not sure how to do a mixture analysis. I guess what I'm asking is how to properly specify the models with the `lcmm()` function. Or, is there maybe another, more convenient way of doing what I want to do in R? Thanks in advance.


title: Error bar for Poisson count data

text: I have a set of data, counts versus time. The whole data looks like this ![data](http://i.imgur.com/XpBbkab.png) Here is a sample                      x     y         [ -0.9   0. ]         [ -0.8   4. ]         [ -0.7   5. ]         [ -0.6   1. ]         [ -0.5   9. ]         [ -0.4   8. ]         [ -0.3  10. ]         [ -0.2   9. ]         [ -0.1  18. ]         [  0.   25. ]         [  0.1  24. ]         [  0.2  26. ]         [  0.3  18. ]         [  0.4  32. ]         [  0.5  25. ]         [  0.6  32. ]         [  0.7  20. ]         [  0.8  23. ]         [  0.9  16. ]      I want to generate suitable error bars for the y-axis in this data. Some people in my research area use the following: $$y \pm \sqrt{y}$$ or $$y \pm y\sqrt{\frac{1}{y} + \frac{1}{\bar{y}}}$$ Because the data is meant to be Poisson, however this doesn't make sense to me because surely for larger y-values, aka more counts...the error should get smaller, not bigger.


title: ADF test in R yielding perfect cointegration. How is this possible?

text: I am using the famous conintegrated pairs tutorial to just different stocks for cointegration. The adf.test yeilds perfect cointegration, which I feel must be incorrect. Here is why: When I run adf.test() on a cumsum of a random series, the plot looks like this:![enter image description here](http://i.stack.imgur.com/KqWda.png) And it yields the following adf.test output: Augmented Dickey-Fuller Test data: sp Dickey-Fuller = -2.8333, Lag order = 4, **p-value = 0.2314** alternative hypothesis: stationary Here is a spread I constructed, notice how it looks similar to the random walk: ![enter image description here](http://i.stack.imgur.com/Q9FZs.png) Which yields the following adf.test() output: Augmented Dickey-Fuller Test data: sprd3 Dickey-Fuller = 3.719, Lag order = 7, **p-value = 0.99** alternative hypothesis: stationary Warning message: In adf.test(sprd3) : p-value greater than printed p-value Any ideas what could be going on here? Why is the p-value extremely different between the two cases? I have a hard time believing that the spread I constructed in the graph is perfectly stationary... Thanks.


title: Interpretation of   "one" feature change in a supervised classifier

text: i'm making experiments using app. 5000 labeled dataset.i'm trying different supervised ML algorithm to evaluate the results.The vector size is 13 with the labels (totally 12 features+1 label) and i have 15 vector of labeled "flower" class. experiments consist of all data set using 10k cross validation. All features are continuous.               1 experiments using the "pure" features of all dataset.     2 experiments using only "one" feature (out of 12) change of the flower class.      i applied naive bayes, C4.5 but all results of 1 and 2 is same, however logistic regression gives different results and lasted longer. 1- To your best experiences, what causes the difference between naive bayes, c.45 and logistic regression, how should i evaluate the results to make the audiences satisfied? 2- if performance is an important metric, and classifier is used for IDS systems, which ML algorithms do you offer? **Edit: More explanation to make the question clear:** we have 8 different class labels. flower + other 7 labels. In experiment 2 we change the only one attribute of the flower class out of 12 attributes (15 flower labeled class stays same but only one attribute is changed.). all other dataset stays same. so we make experiments using logistic regression, naive bayes and c4.5 seperately with two different dataset. (1- with 5000 dataset, 2- other dataset has difference of change in one attribute of flower class, all other classes stays same). comparision: we have results of situation number 1 and 2, in C4.5 and naive bayes, nothing changes, FP and FN. but logistic regression gives interesting results.               12 0.4 0.4 0.5 2.333 434 12.2 10 2 10 12 12 flower      ........................................... flower     ...........................................     ........................................... flower  (total 15 flower class.)          // one feature change 2nd feature.          12 0.8 0.4 0.5 2.333 434 12.2 10 2 10 12 12 flower      ........................................... flower     ...........................................     ........................................... flower  (total 15 flower class,2nd feature all changed.)      > For example can i make i comment like that: because C4.5 uses the maximum > number of class in leaf nodes, the change of one feature in flower class > will not affect/change the leaf node classes However, logistic regression > uses ... so we observe this kind of differences.???


title: Multilevel data with correlation within subjects

text: I have a dataset with 1206 deputies from two different chambers (1998 and 2002, respectively). In addition, there are 18 parties, and some deputies are in both chambers (the ones who were reelected). I am interested in the relation between party discipline (in a scale from 0 to 100) and the probability of geting appointed as party leader (0 or 1), but I am also interested in separate estimates for each party. A classmate told me that the following mixed effect model would be appropriate, but I don't know what to do with the correlation within subjects:               model <- glmer(leader ~ discipline + (1 + discipline|chamber:party),                     data, family=binomial)      What is the structure of these data and which model would be appropriate in this case (using R)?


title: Analysis of temporal patterns

text: I am analysing data on events that I have categorized into groups. So, for instance, say I have 3000 events categorized into 5 groups, which we call A to E. I will have something like                   Event  | Group | Time        --------+-------+-------            1   |   A   |   0            2   |   A   |   5            3   |   C   |   7            4   |   D   |   16           ...  |       |          3000  |   B   |   6000      Now, I would like to see whether there is some sort of n-event long temporal sequence appearing repeatedly (higher than chance). So for instance one 4-event long pattern may be:   **A** \- _3 seconds_ \- **A** \- _2 seconds_ \- **D** \- _5 seconds_ \- **C** I found this paper which proposes some interesting method (I am not actually working on spike trains, but the issue is similar enough), but before implementing that I would like to see whether anyone knew of other methods/statistics that can be applied to this kind of problems.


title: How to estimate errors on a sample with very few data points

text: I have a very simple question to ask, but I can't figure this on my own. I have two samples: sample A with only three data points, and sample B with hundreds of points. For each sample I measure the median of certain quantity. Now my question is, what is the error associated to the median? I can consider the quartiles but, what if sample A consisted of only 1 single point? This would virtually assign no error to median value. I expect my measurement more accurate for sample B.


title: How to calculate disparity of two images in matlab?

text: I have two images and I am trying to calculate the disparity between them using sum of squared distances and reconstruct disparity in 3D space. Do I need to rectify the image before calculating disparity? The following are the steps that I have done so far for disparity map computation(I have tried with rectification and without rectification but both are returning all zeroes disparity matrix).               For each pixel in the left image X,         Take the pixels in the same row in the right image.        Separate the row in right image to windows.        For each window,          Calculate the disparity for each pixel in that window with X          Select the pixel in the window which gives minimum SSD with X        Find the pixel with minimum disparity among all windows as the best match to X      Am I doing it correctly? How can I visualise the 3D reconstruction of the disparity as scatter plot in matlab?


title: Estimating original series from their moving average

text: I have a monthly series of a moving average , and I want to reconstruct their original series from a simple moving average that I have. I happen to know the window of averaging. Am I able to either reconstruct or estimate the original series from their moving average? And in the case I don't have all the info to to fully reconstruct it, how about to estimate it?


title: Non-reverting tendency?

text: I have a time series from which i select a random point. I randomly pick a direction "Up" or "Down" and i place two targets. When i reach either of these targets, i increment the counter for "Target 1 Reached" or "Target 2 Reached" depending on which is reached first. Target 1 is always set at 1 +/- value at the randomly selected point. Target 2 has varying value: 1,2,4,6,8,10 ... +/- value at the randomly selected point. eg. assuming Target 2 = 4               Randomly selected point = 20     Direction = "Up"     Target 1 = 20 + 1 = 21     Target 2 = 20 - 4 = 16      or               Randomly selected point = 18     Direction = "Down"     Target 1 = 18 - 1 = 17     Target 2 = 18 + 4 = 22      If there was a tendency for the time series to remain localized around previous values, shouldnt the "Target 2 Reached" count decrease by more than half, if i double the the size of target 2? Are there betters ways to investigate whether a time series remains localized around certain values?


title: Uniform distribution on 255 from text

text: I'm trying to create a way to link letters from a text to a position between 1 to 255. For example, the text is : "stackexchange" I would like to link every letter to a number between 1 and 255. The same word should have the same output and the overall should be uniform (if I give lots of input, I would like to have an uniform distribution.) // And of course, I want that "stackexchange" always give the same output. // e.g. s=>12 || t=>88 || a=>65 || c=>214 || ... || c=>142 ||... Could you help me ? Thanks !


title: SVM classification step on embedded system with RBF kernel

text: Since it has been some time when i looked into support vector machines and i am facing now to implement the classification step of a trained svm model i would like to ask whether someone could point out how the actual classification step is carried out (assuming i would like to port that step to some low level language) From my trained Matlab SVM model i have: Support vectors (n * #features) the bias (1x1) alpha (n * 1) shift (1 x n) scaleFactor (1 x n) sigma for rbf (1x1) Given a new sample (1 x #features) i would carry out the classification step as follows: Scale and shift each feature in sample: sample = scaleFactor * (sample + shift) calculate the kernel mapping with an RBF with kernel = exp(-1/(2*sigma^2) * ||x-y_i||^2) where x is my new sample and y every single support vector (?) Now i am puzzled, every distance between x and y_i is multiplied by the appropriate alpha? All these values are summed and then the bias is added followed by a simple sign()? So: sign(sum(exp(-1/(2*sigma^2) * ||x-y_i||^2) * alpha_i) + bias) Would that be correct? If so, to save memory on runtime - is there a way to divide the kernel computation in a way that not all support vectors have to be stored in memory? Thanks for your time folks :) Greetings


title: How to test equality of variances with circular data

text: I am interested in comparing the amount of variability within 8 different samples (each from a different population). I am aware that this can be done by several methods with ratio data: F-test equality of variance, Levene test, etc. However, my data is circular/directional... I have done some research and found one test in the "CircStats" package in R - "Watson's test for homogeneity". One shortcoming is that this test only compares two samples, which means I would have to do multiple comparisons on my 8 samples (and then use the Bonferonni correction). Here are my questions: 1) Is there a better test that I can use?   2) If not, what are the assumptions of Watson's test? Is it parametric/non- parametric?   3) What is the algorithm by which I can perform this test? My data is in Matlab, and I would prefer to not have to transfer it into R to run my test. I'd rather just write my own function.


title: Comparison of k-means clustering output

text: I am trying to compare the outputs of k-means algorithm coded by me and the outputs of R's **`kmeans`**. Since the objective of the algorithm is to minimize the total within cluster sum of squares (WCSS), I have to look at the `withinss` output of **`kmeans`** and my WCSS. What I'm not able to figure out is, how much of variation should I allow i.e. how close should these values be so I can say that they are close.   (I know that I can also verify using the centroids, but when the data is not that well clustered, the centroids can be very different but the total WCSS can be close, and hence I'm looking at WCSS.)


title: Showing that E(E(Y|A))|B) = E(Y|A,B)?

text: Would $E( E(Y|A) | B) = E(Y|A,B)$? How does one show this? Intuitively, it make sense that conditioning on this one at a time would result in the same thing as conditioning all at once.


title: Function domain as a problem for linear output unit

text: I'm doing some regression using neural net(using MLP implementation from http://deeplearning.net/tutorial/mlp.html, I used my own but it produced the same results before I opted for this one), standard SGD with backpropagation, learning rate 0.01, weight decay 0.01, activation function for hidden units - $\tanh$, output unit is linear, 10 000 training examples. The function I'm optimising is just $\cos{(x)}$. So the problem is with the domain for x. If it's $[-0.5,0.5]$ then the network just learns straight line: ![enter image description here](http://i.stack.imgur.com/9qT2R.png) (ignore this "predictors" nonsense, lines apart from "nonfunctional" should resemble the target function). If I change nothing but the domain to $[-2, 2]$, we get this: ![enter image description here. -no](http://i.stack.imgur.com/UJH13.png) Which is what we wanted. Why is this happening? I presume that due to some effect of $tanh$ units their input weights quickly become small enough to render all activations completely linear(given that $tanh$ has this linear segment around zero). But why? I seem to be missing something.


title: Relevancy of order statistics to the roll-and-keep dice mechanic?

text: Throw $x$ dice, each of which has $z$ sides, keep the $y$ highest values rolled, and find their sum $s$. In roleplaying games, this is called a "roll and keep dice mechanic" and notated "$x$d$z$k$y$" or perhaps just "$x$k$y$" if the number of sides $z$ is understood. Let (capital) $S$ be the random variable for the sum of the $y$ highest dice, and let (lower case) $s$ be the observed sum of the $y$ highest dice. A combinatorial formula for the probability mass funtion $f$ for this distribution was derived by a user called "techmologist" over at PhysicsForums: Puzzling "roll X dice, choose Y highest" problem. The pmf is $ f(s) = $ $$\sum_{r=1}^{z}\sum_{i=0}^{x-y}\sum_{j=x-y-i+1}^{x-i}{x \choose i,j,x-i-j}(r-1)^i N(s+i+j-x-ry; x-i-j, z-r-1) $$ where $N$ is the number of ways to distribute $B$ indistinguishable balls among $C$ cells so that each cell has between 0 and $D$ balls, namely, $$ N(B;C,D) = \sum_{k=0}^{\lfloor B/(D+1) \rfloor} (-1)^k \binom{C}{k} \binom{B-k(D+1) + C-1}{C-1} $$ My intuition is that there might be a simpler way to express this mechanic using order statistics. Let $Z_1,Z_2,\ldots,Z_x$ be the unordered random variables for the $x$ iid dice throws, which are each discrete uniform distributions on the integer values $\\{1, 2, \ldots, z \\}$. And let $1 \leq Z_{1:x} \leq Z_{2:x} \leq \ldots \leq Z_{x:x} \leq z$ be the order statistics, so that for example, $Z_{x:x}$ is the r.v. representing maximum roll of the $x$ dice. In support of my intuition, it can be shown that $ E(S) = E(Z_{(x-y+1):x} + Z_{(x-y+2):x} + \ldots + Z_{(x-1):x} + Z_{x:x}) $ For instance, for the 4d6k3 distribution, the mean is 15869/1296 (or 12.2446). And this can be verified (quite easily in Mathematica) by using either techmologist's cominatorial formula or by using the expected value of the sum of the three highest order statistics shown above or by enumerating all outcomes or by simulating a large number of dice throws. However, $ S \not= Z_{(x-y+1):x} + Z_{(x-y+2):x} + \ldots + Z_{(x-1):x} + Z_{x:x} $ (WRONG) **EDIT:** Thanks to Douglas Zare for explaining that the LHS and RHS are in fact equal... This helped me realize that my real confusion was that the pmf of S cannot be calculated easily by taking the convolution of the pmfs of the second through fourth order statistics, because the order statistics are not independent RVs. And that leads to a follow-up question: In a case like this, how would you find the convolution of the pmfs of the order statistics on the RHS, given that order statistics are not independent RVs? I will think about a better way to phrase this question and post it separately. Thanks again! The question is whether my intuition is correct: **Is there some way to express the pmf $f$ of $S$ in terms of the pmfs of the order statistics? And if not, what's wrong with my intuition that there should be some easy way to calculate the pmf $f$ from the pmf of the order statistics?** BTW, this isn't homework. I'm just a roleplayer, not a statistician, trying to satisfy my curiosity whether this dice mechanic can be analyzed and expressed more simply using order statistics. Thanks.


this is a somewhat random question: is it possible to find clustered standard errors for first difference regressions? (both theoretically and practically in R or Stata). Thanks!


title: Binary choice model and classification

text: Are there any differences between binary choice model and classification problem? For example, binary choice model may use logistic regression, which is also an widely used approach in classification.


title: Iterating through rows and filtering data in R

text: I am working with the following data (see picture) I am trying to predict future Y axis values from the X values. As you can tell from the picture and as I discovered, the prediction is going to be horrible based on the data alone. I am now trying to go row by row through the data and classify each point as significant or nonsignificant and see if I can find anything interesting in the classification. I want my algorithm something like this:   1. Take absolute value of every element in both Y axis and X axis data columns   2. Iterate through Y axis column and find the top 5% largest values   3. Iterate through again marking those top 5% or maybe 1% (Will fiddle with this) with a flag   4. Plot only the top 5% and see if the linear model (or some model) has a better R^2   5. Perform some sort of classification analysis against the two flags, significant and non significant (I really need help and suggestions with this step) ![enter image description here](http://i.stack.imgur.com/renrN.png) The two things I need are a pointer on are   * How to step element by element through data, making modifications along the way   * What kind of classification analysis exists that might be helpful, goal being to predict which classification future sampled elements will fall into. Thanks!


title: What statistical test do I need to see if there is a correlation between two different instruments?

text: I want to do a study to examine if there is a correlation between nursing students' confidence & satisfaction, and their intuitiveness. I have two instruments, _The Student Satisfaction and Self-Confidence in Learning questionnaire_ is a 13 item instrument which uses a five point scale to measure self confidence and satisfaction in nursing students and the _Miller Intuitiveness questionnaire_ , which is an 18 item instrument which uses a six point scale. What statistical test can I use to see if there's a connection between the two? Edit: I'm leaning towards Spearmans Rank Order correlation, can anyone confirm this is suitable please?


title: How to best display crosstab data?

text: I have a 10x10 matrix composed of two variables with 10 brands each. One variable is the brand purchased, the other is the brand considered. My matrix shows a crosstabulation between the two. I need an effective way to clearly visualize these data so that proximities suggest similarities between categories and distances dissimilarities. I ran a correspondence analysis, but wasn't impressed with the graph. Are there any alternative technique to consider?


title: Validation for ANOVA question

text: Can someone help me validate my answer to this question: ![enter image description here](http://i.stack.imgur.com/CwlLR.png) My answers come from the first row, where the significances are 0.0002,0.0011, and 0.0005, respectively. Since they are all below 5%, then the conclusion should be correct, right? Am I thinking about this incorrectly?


title: Estimate specific y value in linear multiple regression using R

text: I am trying to find specific y value of fitted regression model using specific x value. for example, in the multiple regression model like below, > y~ aX+bZ+cW +error **I fitted the model with data** and I have the result of fitted glm model like > y = 3.142X+2.94Z+4.16W + 2.194 I'm trying to get 'y' when X is 2 and X is 5. Is it possible in R?


title: Derivation of bayesian confidence interval

text: How do you derive the expression for the $100(1-\alpha)$% bayesian confidence interval when working with the uniform distribution in the interval $[-\theta,\theta]$?


title: How to avoid exact fit at entry level in decision tree?

text: I implemented decision tree according by Russell and Norvig AIMA (in the 3rd edition, p.704). ## Assumption The given formula for calculating entropy does not work for the cases when you have all positives examples, or all negative ones. For binary outcomes: $$B(q) = -(q * lg_2(q)+(1-q) * lg_2(1-q))$$ $$q=\frac{p}{p+n}$$ > `p` denotes number of positives, `n` negatives So I made a correction, that in such cases nothing is computed, and simply zero (total order) is returned. ## Exact fit problem Let's say you have info about population in USA, and you would like to filter out all people living in Boston. However if you have in your attributes social number, you could assign "Boston" (or "not Boston") with pin-point accuracy. So when DT computes gain per each attribute, when it comes to social number, all of the sudden gain is maximum, because each attribute value (particular social number) gives immediate outcome about location. Let's take a look at the equation (I rewrite equation for remainder, maximum gain = minimum remainder): $$ Remainder(A) = \sum_{k=1}^d \frac{p_k+n_k}{p+n}B(\frac{p_k}{p_k+n_k}) $$ > B denotes entropy, `k` indicates given class (value of attribute) So, having `n` people DT divided entire population into `n` subsets, each of 1 person. It was possible because entropy of such subset is zero, product is thus also zero, and the sum in effect is zero. I will have `n` edges from root node (an overkill and overfit at the same time), yet for DT it is optimal. ## Question How to prevent such effect? AIMA touches this problem too lightly for me (p.707). ### My guess Every multi-value attribute convert into binary-value set of attributes. For example for colors -- convert:               color = { red, orange, blue }      into               is_red = { true, false }     is_orange = { true, false }     is_blue = { true, false }      This way I could avoid zero-entropy problem, and if I finally get to such situation it would be for good reason. Second, I could be sure about number of edges, so implicitly this factor would be embedded into implementation (in AIMA this is out of control, see the equation above, number of edges -- `k` \-- is not used as a factor). Another approach would be adding `k` as an entropy factor of edges. However I have no clue how to do it reliably.


title: How do to calculate Likelihood Ratio Test/Power in hypothesis testing?

text: I am preparing for an exam, I've been reading through Likelihood Ratio Test, but don't get it. Example 8.2.3 (page 376 in "Statistical Inference" Casella and Berg) Let ${{x}_{1}},...{{x}_{n}}$ be a random sample from an exponential population with pdf $$f(x|\theta )=\left\\{ \begin{array}{*{35}{l}} {{e}^{-(x-\theta )}} & x\ge \theta \\\ 0 & x<\theta \\\ \end{array} \right.$$ ${{H}_{0}}_{{}}:\theta \le \theta $ versus ${{H}_{1}}:\theta >{{\theta }_{0}}$ The Likelihood function is, $$L(\theta |x)=\left\\{ \begin{array}{*{35}{l}} {{e}^{-\Sigma {{x}_{i}}+n\theta }} & \theta \ge {{x}_{(1)}} \\\ 0 & \theta <{{x}_{(1)}} \\\ \end{array} \right.$$ The book says: Clearly, $L(\theta |x)$ is an increasing function of $\theta$ on $-\infty < \theta \le {{x}_{(1)}}$. (Does this even matter?). Thus, the denominator of $\lambda(x)$, the unrestricted maximum of $L(\theta|x)$, is $$L({{x}_{1}}|x)={{e}^{-\Sigma{{x}_{i}}+n{{x}_{(1)}}}}$$ If ${x}_{1}\le{\theta}_{0}$, the numerator of $\lambda(x)$ is also $L({x}_{(1)}|x)$. But since we are maximizing $L(\theta|x)$ over $\theta\le{\theta}_{0}$, the numerator of $\lambda(x)$ is $L({\theta}_{0})$ if ${x}_{1}>{\theta}_{0}$. Therefore, the likelihood ratio test statistic is $$\lambda (x)=\left\\{ \begin{array}{*{35}{l}} 1 & {{x}_{(1)}}\ge {{\theta }_{0}} \\\ {{e}^{-n\left( {{x}_{(1)}}-{{\theta }_{0}} \right)}} & {{x}_{(1)}}<{{\theta }_{0}} \\\ \end{array} \right.$$ Rejection region $$\\{x:{{x}_{(1)}}\ge {{\theta }_{0}}-\frac{\log c}{n}\\}$$ I read through the whole thing, it was long and tedious to follow, I'm quite lost. Can someone else explain this? Maybe with a better/simple toy example? I want to understand this, not just follow through! Is this a simple concept, I am over-complicating it?


title: Comparing two discrete distributions (with small cell counts)

text: I need to compare sample distributions with theoretical ones, which is typically done with a chi-squared test. The problem is that I have distributions where one or more cells have low values, and consequently the chi-squared test reports very small p-values. For example, a typical expected and observed frequencies are [152 2 9] and [140 5 18], with a p-value of 0.0007. Based on domain knowledge, these two distributions are not significantly different. What test could be used instead of chi-squared, which would take out the bias that occurs with the small-valued cells? **Edit: adding some background information for this problem.** I have a number of processes which produce as output certain technical parameters, recorded as time series. I have around 4000 of such process, each producing around 150 such time series (the number of time series a process has follows a power law). I would like to find which of these processes are anomalous, i.e. producing output which is significantly different from others. To do this I cluster the time series using k-means, and then based on the clusters, produce the "expected" distribution (average over all time series) and the distribution of clusters for each process. For example, after clustering I might have 4 clusters with following sizes.               Cluster number | Cluster size     -----------------------------     1              | 100     2              | 200     3              | 300     4              | 400      The distribution of the clusters among the processes might be the following                         | Cluster 1 | Cluster 2 | Cluster 3 | Cluster 4     ----------------------------------------------------------     Process 1 | 11        | 19        | 35        | 42     Process 2 | 3         | 10        | 14        | 19     Process 3 | 30        | 8         | 12        | 12              <----anomaly     ....      In this case, process 1 and process 2 are sufficiently close to the expected, while process 3 has a different distribution from the average. I would like to find a good test to measure this discrepancy. (Any other suggestion for the anomaly detection is also welcome)


title: How to test a curvilinear relationship in a logistic regression

text: I was looking for some information about curvilinear relationships (quadratic function, to be precise) in logistic regression online, and couldn't really find much about it. I am interested if that could be done, and if so, how can I do it in SPSS? Do I just square the variable that I want to test for a curvilinear relationship and add it into the model like any other variable? Or is there some special way of doing it?


title: Is a p-value of 0.04993 enough to reject null hypothesis?

text: In a Wilcoxon signed-ranks statistical significance test, we came across some data that produces a $p$-value of $0.04993$. With a threshold of $p < 0.05$, is this result enough to reject the null hypothesis, or is it safer to say the test was inconclusive, since if we round the p-value to 3 decimal places it becomes $0.050$?


title: transformation of percentages (and eventual display)

text: I have some percentage data derived from the analysis of grain presence and absence in ears of wheat that have undergone two treatments (control and heat stress). For example, data might be like 10%, 15%, 20% for stress and 80%, 85%, 90% for control. I understand (not least from the link below) that I need to transform my data to arcsin and then analyse it as you would any other data set. This is correct, isn't it? http://archive.bio.ed.ac.uk/jdeacon/statistics/tress4.html#Converting My primary question is once I have transformed this data, and stated that I have done so in my methodology (e.g. all other percentages were subjected to angular transformation and these, along with the rest of the data, were analysed by parametric statistical methods), do I continue to refer to the data in graphs and tables as %? For example, does 8%, once transformed, become arcsine 16.48 or can I still describe it as 16.48% and have % on the y axis of a graph?


title: How to find Beta?

text: I have a dataframe (of prices) as given below               dat <- data.frame(company = rep(c("A", "B", "C", "D", "index"), each = 5),     prices = c(runif(5, 10, 12), runif(5, 108, 112), runif(5, 500, 510),     runif(5, 40, 50), runif(5, 1000, 1020)))      I need to find the beta of A, B, C and D w.r.t index. Beta between two variables X and Y (where Y is dependent) is given by               beta <- coef(lm(Y ~ X))[2]      Any guidance is appreciated.


title: how can I calculate $cov(T_1,T_2)$

text: Suppose $X_1,\ldots,X_n$ is a random sample from the density function $$f(x,\theta)=\displaystyle\frac{1}{2\theta}e^{\displaystyle-\frac{|x|}{\theta}}, \quad -\infty<x<\infty\,\text{.}$$ If $T_1=\sum_{i=1}^{n} a_i X_i\, $ and $T_2=\sum_{i=1}^{n} |X_i|$, how would I calculate $cov(T_1,T_2)$?


title: Does transformation of r into Fisher z benefit a meta-analysis?

text: Usually $r$ is transformed into Fisher $z$ to test difference between two $r$ values. But, when a meta-analysis is to be performed, why we should take such a step? Does it correct for measurement error or non-sampling error and why should we assume that $r$ is an imperfect estimate of population correlation?


title: Explain probability measure and probability distribution

text: Dear statistics community, I seek your help to understand the jungle of concepts about probability theory. I am struggling to form a coherent understanding of terms like _probability mass function_ , _probability distribution_ , _probability function_ , concepts about random variables and the ambiguous notation surrounding them. E.g. both german and english Wikipedia avoid to explain the relationship between   1. the function mapping events of a probability space (this is rather consistently defined) to numbers, sometimes calling those numbers probabilities and   2. the concept of probability distribution. Is the function mapping events (called a _probability measure_ by en Wiki) to numbers a _probability distribution_? If not what is a probability distribution?


title: How to compare response to same question asked at different time points?

text: I asked survey respondents about two points in time, and asked participants the degree to which they "agree", a typical Likert-style question. For example, as if I had asked "Before Obama became president, life was good" and a separate question, "After Obama became President, life has been good". Normally, one would assume a t-test, and at first I thought a paired t-test, since these are two separate sets of answers I am trying to compare (sort of like time 1/time 2, right?). However, I've gotten rejected multiple times by wanna-be stats folks, none of whom can tell me what the right test is. So can anyone tell me if I ask about two points in time (in your childhood, what was your relationship like with your mother? Now, what is your relationship with your mother?), how do I compare them?


title: Imputation variance and explained variance (in vector autoregression)

text: I have a question concerning the coefficients of VAR models used on multiple imputed data (high missigness in some variables: up to 40%). In particular I would like to know how the coefficients are related to the explained variance. I have used vector autoregression on multiple imputed data (m=10) and have then combined the estimated coefficient with rubin's rule. However, what confuses me is the fact that my imputation variance is quite small in relationship to the estimates and variance of coefficients, but the difference between the explained variance is huge (17% to 0.04%) between models. My idea is that since the highest imputation variance across all systems is at the constant (around a third of the variance value but 3-4 times higher then in other coefficients) and that this critically affects the explained variance. But thats just a guess. I would be very happy if somebody could help me here.


title: Transformed data odds ratio interpretation in a logit model

text: One of the independent variables for my binomial logit regression is age. In order to make age normal, I transform it to square root age by using ladder options. How I can interpret the transformed square root age using odd ratio. Thank you!


title: Distribution of a ratio of uniforms: What is wrong?

text: Suppose that $X$ and $Y$ are two i.i.d. uniform random variables on the interval $[0,1]$ Let $Z=X/Y$, I am finding the cdf of $Z$, i.e. $ \Pr(Z\leq z) $. Now, I came up with two ways of doing this. One produces a correct answer consistent with the pdf here: http://mathworld.wolfram.com/UniformRatioDistribution.html, the other does not. Why is the second method wrong? **First Method** $\newcommand{\rd}{\mathrm{d}} \Pr(Z\leq z) = \Pr(X/Y\leq z) = \Pr(X\leq zY) = \int^{1}_{0}\int^{\min(1,zy)}_{0} \rd x \rd y = \int^{1}_{0}\min(1,zy)\ \rd y$ $ = \left\\{ \begin{array}{lr} \int^{1/z}_{0}zy\ \rd y + \int^{1}_{1/z} \rd y& : z > 1\\\ \int^{1}_{0}zy\ \rd y & : z \leq 1 \end{array} \right. $ $ = \left\\{ \begin{array}{lr} 1 - \frac{1}{2z} & : z > 1\\\ \frac{z}{2} & : z \leq 1 \end{array} \right. $ This appears correct. **Second Method** $\Pr(X/Y\leq z) = \Pr(X \leq zY\ |\ zY \geq 1)\Pr(zY \geq 1) + \Pr(X \leq zY\ |\ zY < 1)\Pr(zY < 1) $ by total probability $ = \Pr(X \leq zY\ |\ zY \geq 1)\Pr(Y \geq 1/z) + \Pr(X \leq zY\ |\ zY < 1)\Pr(Y < 1/z)$ Taking $z>1$ yields $ (1)(1-\frac{1}{z}) + (\int^{1/z}_{0}\int^{zy}_{0} \rd x \rd y)(\frac{1}{z}) = 1-\frac{1}{z} + (\int^{1/z}_{0}zy\ \rd y)(\frac{1}{z}) = 1-\frac{1}{z} + \frac{1}{2z^{2}}$ This is already different. Why is this wrong? Thanks!


title: Help with Euclidean update rule

text: I'm trying to figure out the derivation of the Euclidean update rule for NMF on page 2 of http://sig.umd.edu/publications/Tjoa_ICASSP2_201003.pdf but can't seem to eke it out. The end result that they get is: $$ \frac{\partial}{\partial \textbf{S}} d_{EUC}(\textbf{Q}, \textbf{S}\textbf{S}^{T}) \propto \textbf{S}\textbf{S}^{T}\textbf{S} - \textbf{Q}\textbf{S} $$ my working was to let $\textbf{X} = \textbf{Q} - \textbf{S}\textbf{S}^{T}$ and thus $$ \frac{\partial \textbf{X}}{\partial \textbf{S}} = -(\textbf{S}^T + \textbf{S}) $$ so using the chain rule $$ \frac{\partial}{\partial \textbf{S}} d_{EUC}(\textbf{Q}, \textbf{S}\textbf{S}^{T}) = \frac{\partial}{\partial \textbf{S}} \text{tr}((\textbf{Q} - \textbf{S}\textbf{S}^T)^T(\textbf{Q} - \textbf{S}\textbf{S}^T)) = \frac{\partial}{\partial \textbf{X}} \cdot \frac{\partial \textbf{X}}{\partial \textbf{S}} \text{tr}(\textbf{X}^T \textbf{X}) $$ and since $\frac{\partial}{\partial \textbf{X}}\text{tr}(\textbf{X}^T \textbf{X}) = 2\textbf{X} $ this reduces to: $$ \frac{\partial}{\partial \textbf{S}} d_{EUC}(\textbf{Q}, \textbf{S}\textbf{S}^{T}) \propto (\textbf{S}^T + \textbf{S})(\textbf{S}\textbf{S}^T-\textbf{Q}) $$ which is clearly wrong. I'm not too sure where it all goes awry, I think there's probably some very debatable maths in there. If anyone can steer me in the right direction I'd be very grateful!


title: Comparing regression outputs for different response variables

text: I have inherited some legacy code designed to solve the following problem - we are given a set of observable binary variables $Y_1, \dots , Y_N$ that we believe are related to a fixed set of input variables $X_1 ,\dots , X_M$ (N.B. high dimension, $M \sim 100$). For a given input we want to predict which event $Y_i$ is the most likely. The current program tries to predict this as follows - given training data, we build a model for each $Y_i$ seperately (using logistic regression + feature selection). Then given an input we look at all of the logit scores and choose the highest. This doesn't work because some of the logistic regression classifiers can be overtrained, which leads to extreme logit scores which often dominate the other scores. Analysis reveals that the logistic function may not be a good choice either. I have two questions :   * is this way of finding the most likely event sensible / salvagable? I would rather not have to rebuild this whole thing from scratch if possible.    * if so, is there a more sensible way to compare the scores from two classifiers in such a way that it is more fair? My thoughts are that we could compare some sort of Bayesian posteriori to get a fair evaluation for our state of knowledge, but I feel like if I were to do this I would have to lose the logistic regression model, and I am not sure what I should replace it with.


title: Mann whitney U test

text: I would like to compare the nutrient intake of males with that of females. One of the assumptions of the Mann Whitney test U test is that the two variables must not be normally distributed. The nutrient intake of males is normal but not of females. Would this be a problem? Also, the data is on intake of seven nutrients. To do the comparison wouldn't I have to do the Mann whitney test 7 times? Is there any other alternative? The following are the qq plots I obtained for women and men respectively. ![Chi square qq plot for Women nutrient intake](http://i.stack.imgur.com/Sxjhh.jpg) ![Chi square qq plot for men nutrient intake](http://i.stack.imgur.com/VhWfp.jpg)


title: Good Adjusted r square

text: what is a good adjusted R square value for any multiple regression model or it depends? If it depends then what are the conditions?


title: Bias Variance Decomposition for Mean Absolute Error

text: The mean squared error of an estimator $\hat{\theta}$ with respect to an unknown parameter $\theta$ is defined as $$ MSE(\hat{\theta})=E[(\hat{\theta}-\theta)^{2}]. $$ It is well known that there is the following bias-variance decomposition: $$ MSE(\hat{\theta})=Var(\hat{\theta})+\left(Bias(\hat{\theta},\theta)\right)^{2}. $$ **My question** : does the mean absolute deviation error $$ MADE(\hat{\theta})=E(|\hat{\theta}-\theta|) $$ have a similar decomposition?


title: How Should I Generate a Plausible PDF From Quartiles+Median?

text: I'm developing a tool to conduct expert elicitations of scientists. As answers to questions are speculative, we like experts to help us understand the range of answers they find plausible, to get some idea of uncertainty. For one type of question, we have them input a lower-bound, upper-bound, 25th percentile, 50th percentile and median (essentially the "5 number summary"). In testing we found that some scientists find this data intuitive, and others find it easier to visualize in the form of a PDF. I'm trying to figure out a straightforward algorithm to provide a PDF-style visualization of this "5 number summary". I've considered a "Skew Normal" distribution, would this be a good choice? Put another way: given a quartiles+median, what's the most straightforward function to render a PDF for a distribution with those characteristics? Any suggestions? I know this is sloppy at best, but the answers to these questions are speculative themselves, and we'll get better answers if we can communicate both ways.


title: Is there a theorem that says that $\sqrt{n}\frac{\bar{X} - \mu}{S}$ converges in distribution to a normal as $n$ goes to infinity?

text: Let $X$ be any distribution with defined mean, $\mu$, and standard deviation, $\sigma$. The central limit theorem says that $$ \sqrt{n}\frac{\bar{X} - \mu}{\sigma} $$ converges in distribution to a standard normal distribution. If we replace $\sigma$ by the sample standard deviation $S$, is there a theorem stating that $$ \sqrt{n}\frac{\bar{X} - \mu}{S} $$ converges in distribution to a t-distribution? Since for large $n$ a t-distribution approaches a normal, the theorem, if it exists, may state that the limit is a standard normal distribution. Hence, it would seem to me that t-distributions are not very useful - that they are useful only when $X$ is roughly a normal. Is this the case? If it is possible, would you indicate references that containing a proof of this CLT when $\sigma$ is replace by $S$? Such a reference could preferably use measure theory concepts. But anything would be great to me at this point.


title: Generating data with given regression coefficients

text: I would like to generate data points with given regression coefficients, where some of the variables may be discrete and some continuous, and I am performing both "direct" and "reverse" regressions. For example, I would like to generate data points consisting of triples $(x,y,z)$ where $x$ is either 0 or 1 and $y$ and $z$ are real, such that when one (linearly) regresses $z$ on the other two variables ($z_i=\alpha x_i+\beta y_i+\epsilon_i$), I have some given regression coefficients ($\alpha$ and $\beta$). Simultaneously, I would like to specify the coefficients when one regresses y on the other two variables. (As an added complication, I would also like to specify the regression coefficients of $z$ and $y$ conditioned on the different values of $x$.) The difficulty arises, naturally, due to the fact that I would like control of both of these regressions simultaneously. How would I go about doing this systematically?


title: Reporting significance in Kruskal Wallis test

text: I performed the Kruskal-Wallis one-way ANOVA by ranks on five groups of data using Dataplot software. The software corrects for ties and uses Conover's procedure to do pairwise multiple comparisons. The result of the test is $H=5.62$ and $p=.229$ so the results are not significant ($p$ is not less than $.05$). However, the pairwise multiple comparisons show that the group 1 data significantly differs ($p<.05$) from the group 2 data. I report that there are no significant differences in the five groups (because $p=.229$). My question is this: I am using box plots to graphically represent the results. Should I indicate that the group 1 data and the group 2 data significantly differ even though the H-test says that there are no significant differences in the five groups of data?


title: Identifying outliers

text: I have binomial frequency data for an allele associated with populations living in mountainous. These mountains run north to south where sites are nearly fixed for this allele and lowland sites to the east that lack the allele. At the south end of the mountains is a hybrid zone. Despite the mountains continuing south through the hybrid zone, the allele frequency shifts from nearly fixed to nearly absent. This suggests that despite the assumed benefit for this trait in the mountains, that there is a barrier for it in the hybrid zone preventing its spread south in the mountains. I have used a simple glm approach in R to identify that altitude is a significant predictor of the allele frequency north of the hybrid zone. If I include hybrid zone sites in the model, altitude is less significant. If I add a hybrid zone dummy variable to the model and it is also significant. I would like to find a statistical approach to identify these "outlier" mountain sites in the hybrid zone both statistically and graphically.


title: Is it meaningful to test equality of two coefficients, one significant, the other insignificant?

text: My question is like this. The regression specification is: $y=\delta_{1}D_{1}+\delta_{2}D_{2}+\delta_{3}D_{3}+\beta_{1}x_{1}+\beta_{2}x_{2}+\beta_{3}x_{3}+\epsilon$ where $D_{i}$ denotes dummy for each one of three categories, $x_{i}$ interaction between $D_{i}$ and independent variable $x$. I want to test $b_{1}=b_{2}=b_{3}$. $b_{i}$ is the estimated coefficient of $x_{i}$. After test, I will combine interaction terms for the categories not showing significantly different coefficients. Each category has different sample size. Category 1 have much more samples than category 2, and category 2 more than category 3. Estimates show that $b_{1}$ and $b_{2}$ is significant (from zero), but $b_{3}$ not. I think the small sample size of category 3 could explain insignificance of $b_{3}$. My question: Is it meaningful (or right) to do joint test $b_{1}=b_{2}=b_{3}$ when $b_{3}$ is insignificant? If the test $b_{1}=b_{2}=b_{3}$ can not be rejected at e.g 5% level but $b_{1}$ and $b_{3}$ is significantly different, could I say no significant difference among $x_{1}$, $x_{2}$ and $x_{3}$?


title: Fisher information $J_y(\theta)$ for transformation $y=F(x)$

text: Consider a multivariate random variable $x$ with density function $P_x(\theta)$ for a scalar parameter $\theta$. Assume the Fisher information $J_x(\theta)$ is known. Now, for a transformation (deterministic mapping) $y=F(x)$ can you express $J_y(\theta)$ in terms of $J_x(\theta)$ and $F$?


title: How to incorporate sample groups?

text: The regression specification is: $y=\delta_{1}D_{1}+\delta_{2}D_{2}+\delta_{3}D_{3}+\beta_{1}D_{1} \cdot X+ \beta_{2}D_{2} \cdot X+\beta_{3}D_{3} \cdot X+\epsilon$ where $D_{i}$ denotes dummy for each one of three groups, $D_{i} \cdot X$ interaction between $D_{i}$ and independent variable $X$. Each group has different sample size. Group 1 has much more samples than group 2, and group 2 more than group 3. If I can incorporate some of these groups, I will remedy the problem of small sample size in group 2 and 3 to some extent. So I want to test $\hat{\beta_{1}}=\hat{\beta_{2}}=\hat{\beta_{3}}$. $\hat{\beta_{i}}$ is the estimated coefficient of $D_{i} \cdot X$. After test, I will combine interaction terms for the groups not showing significantly different coefficients. The questions are:   1. Should I use joint test $\hat{\beta_{1}}=\hat{\beta_{2}}=\hat{\beta_{3}}$? If this test cannot be rejected at e.g. 5% level, I will incorporate group 1, 2 and 3. Is that right?    2. If I should use joint test, since group 2 and 3 has very few samples, $\hat{\beta_{2}}$ and $\hat{\beta_{3}}$ may be insignificant even having wrong sign (I mean they are should be positive, but estimates are negative.) Is it right to use joint test $\hat{\beta_{1}}=\hat{\beta_{2}}=\hat{\beta_{3}}$ here?   3. Or should I use pairwise tests? If I use pairwise tests to find that $\hat{\beta_{1}}=\hat{\beta_{2}}$ and $\hat{\beta_{2}}=\hat{\beta_{3}}$ cannot be rejected but $\hat{\beta_{1}}=\hat{\beta_{3}}$ is rejected, how to incorporate these groups?


title: Is there overfitting in my modelling approach despite cross-vaidation?

text: My model is predicting a binomial dependent variable with a rich feature space of 20,000 independent variables. I am using the penalized logistic regression from the glmnet package, which works for the kind of dataset I have. I am using 10-fold cross-validation and the cross validation error is just 2.4%. I was told that with such a small cross validation error (in other words a very high accuracy of 97.6%) there is a high chance that my model is overfitting. Is is true? If yes, then what is an alternative approach to this? Here is some code which I am using:               library('glmnet')     data <- read.csv('datafile.csv', header=T)     mat = as.matrix(data)     X = mat[,1:ncol(mat)-1]      y = mat[,ncol(mat)]     fit <- cv.glmnet(X,y, family="binomial", type.measure = "class")     betacoeff = as.matrix(fit$glmnet.fit$beta[,ncol(fit$glmnet.fit$beta)])


title: Why is 2-fold CV giving better results than 5-fold CV?

text: I am working on a binary classification problem with relatively few instances (e.g. ~30 instances out of which ~7 are positives). I have noticed that when using 2-fold the average classification performance of the best performing model is better than the best performing model with 5-fold. In fact,   * The best performing model in 2-fold CV gets the following scores across the two folds: `[0.82, 0.82]` ( **avg. = 0.82** ).   * That model is **different** from the best one I get with 5-fold CV, which yields the following AUC scores: `[0.4 , 1. , 0.75, 0.75, 0.25]` ( **avg = 0.64** ). This takes me to the following question: Which model should I use? And why would I ever get a better model when training with **less** data?


title: Test for aggregation of binary events/successes (binomial/glm??)

text: this has been vexing me for a while and I can't seem to solidify an answer beyond vague thoughts about Poisson distributions. I think this is a simple problem and I'm missing something obvious. Any thoughts appreciated, particularly with R snippets :) Chris Experiment: 2 replicate assays from each of 50 tissues. In the 100 samples, call a binary (0,1) event (without using tissue of origin; I'll call 1 outcomes "successes", per Bernoulli trial convention). Then score tissues for concordance (2 if both replicates are successful; 1 if only one; 0 otherwise). If the assay is working, we expect non-random numbers of concordance = 2 for tissues (i.e. the effect should replicate), given the total number of successes S. How can I test that the distribution of (0,1,2) on _tissues_ \- or perhaps the proportions of (1,2) _tissues_ \- deviates from expectation? thanks!


title: R AUC never less than 0.5?

text: I'm doing some work with random forests in R using the `randomForest` package, and I've run into something that seems odd to me. Even when the data is completely random, the AUC is never less than 0.5. For example, when I run the following:               library(randomForest)     df.sanity <- data.frame(A=sample(1:100, 2000, replace=T), B=sample(126:159, 2000, replace=T), C=sample(10:2000, 1000, replace=T), D=sample(1:2, 2000, replace=T), E=sample(30:40, 2000, replace=T), Class=as.factor(sample(0:1, 2000, replace=T)))     rf <- randomForest(x=df.sanity[1:1000,c("A", "B", "C", "D", "E")], y=df.sanity[1:1000, "Class"])     preds <- predict(rf, newdata=df.sanity[1001:2000,], type="prob")     auc(obs=df.sanity[1001:2000, "Class"], pred=preds[,2])      No matter how many times I run it, the AUC is never less than 0.5. It's often a bit over (up to 0.54 from what I've seen), but never less. The only other AUC implementation I've used is Weka's, and I've seen AUCs < 0.5 there. Does the `randomForest` package automatically flip the predictions to the reverse if the AUC is ever less than 0.5, or is there something else I'm misunderstanding here?


title: Ordinal Logistic Regression with a Different Link Function

text: Consider an outcome variable that has four clear, ordered categories to it. This seems like a good use of ordinal logistic regression to estimate Odds Ratios for the effect of covariates on moving a subject one "step" up the ladder. But the subjects are particularly evenly spread throughout the categories, so a question arises:   * Is the "rare outcome assumption" for an OR to approximate a relative risk still true in ordinal logistic regression?   * If so, is it possible to change the link function to directly estimate a relative risk, and is it still possible to use something like a poisson approximation with robust standard errors to deal with convergence issues in such a case?


title: Stratified random sampling and its distribution

text: I'm learning stratified random sampling and I'm having hard time to solve this problem. consider the following population of six units:               x1=1, x2=2, x3=2, x4=4 x5=4 x6=5      Suppose units 2 3 4 and 5 are in one stratum and units 1 and 6 are in a second stratum. Take a simple random sample of 2 units from the first stratum and a simple random sample of 1 unit from the second stratum. What is the exact distribution of the stratified estimator for the population average? So, possible mean of sample from first stratum is 2, 3 or 4 and from second stratum is 1 or 5. overall sample mean is               2/3 * 2 + 1/3 * 1 = 5/3     2/2 * 2 + 1/3 * 1 = 7/3     2/3 * 4 + 1/3 * 1 = 3     2/3 * 2 + 1/3 * 5 = 3     2/3 * 3 + 1/3 * 5 = 11/3     2/3 * 4 + 1/3 * 5 = 13/3      and each probability is `1/12, 4/12, 1/12, 1/12, 4/12, 1/12`. Now what should I do for finding exact distribution?


title: Resources for learning statistics - exercises (with solutions), available online?

text: I'm currently working as a teaching assistant at my University, in an introductory statistics course (for medical students). Offline, there are many books available with information to aid the teacher. However what I am interested to know is if you might direct me to any (good) **resources that provide exercises (with solutions) in statistics, which are available online?** (e.g: teachers notes). The subject material can range between descriptive statistics, probability, and parametric/non-parametric statistical inference.


title: When to use CDF and PDF for Exponential Distribution

text: I am doing a self-study on Exponential Distribution and noted an exponential on my text giving that: CDF of Exponential Distribution $$ F(x) = 1 - e^{-λx} , $$ PDF of Exponential Distribution $$ f(x) = λe^{-(λx)} . $$ However, I was wondering on what conditions do I use what? From what I understand, if I was trying to find the time between consecutive events within a certain period of time, I may use the CDF. However, I am unable about PDF. Appreciate any advice please.


title: Is there an unpaired version of the sign test?

text: I wasn't able to find one on Google, and was wondering if one exists...


title: What are the sharpest known tail bounds for $\chi_k^2$ distributed variables?

text: Let $X \sim \chi^2_k$ be a chi-squared distributed random variable with $k$ degrees of freedom. What are the sharpest known bounds for the following probabilities $$ \mathbb{P}[X > t] \leq 1 - \delta_1(t, k) $$ and $$ \mathbb{P}[X < z] \leq 1 - \delta_2(z, k) $$ where $\delta_1$ and $\delta_2$ are some functions. Pointers to relevant papers would be appreciated.