title: Is this flawed statistics?

text: I often see similar statements on football (soccer) websites. > Team A has lost 4 out of their last 5 matches against Team B. (and thus, > there is a high chance Team A will lose their next match tomorrow against > Team B. But if I procede to find out the matches played, the most recent of their last 5 matches would be something like 10 years ago! So in that span of 12 years (between 2000 and 2012), the entire team would have changed, the manager has changed etc. So is it flawed to write statements like that? * * * Also, if a team just plays 1 game (assuming you can only win or lose, no draw), looking at their win percentage without seeing the number of games played could be misleading as the win percentage could be 100% or 0%. But as they play more games, you get a better idea of their "actual" win ratio. Is there some sort of statistical term for that?


title: The test procedure of significant differences in count data over time

text: The data below are from a health testing centre over 8 years so there are different people being tested there each year. All I really want to show is that the amount of females being tested is different over the different years. The data look like this:               YEAR   FEMALE   MALE   TOTAL     2002     50       7      57     2003    126      50     176     2004    339      68     407     2005    401     158     559     2006    639     328     967     2007    744     600    1344     2008    570     526    1096     2009    452     448     900     2010    561     372     933      I think these data are not actually repeated measures, because each year different individuals are being tested. Moreover, the data are counts in nature. So, I was just thinking of using Runs test to show if the amount of females being tested is different over the different years. Because Runs test is a test of randomness. The alternative hypothesis should be indicating some trend in the data. Am I right that the data are not repeated measures? Moreover, is the concept of Runs test appropriate here? What would the test be if I wanted to test "if the proportion of females being tested is different over the different years"? I think this is a more sensible test to do, but how can I do it? I have several comparisons like this, so it would be good to know if I'm using the wrong test! Any suggestion regarding the test procedure will be greatly appreciated.


title: How is Item Response Theory (IRT) called for continuous response?

text: I would like to model my problem using something similar to Item Response Theory - but my responses are not binary - they are continuous in [0;1]. How are these models/the research field called?


title: The relationship between expectation-maximization and majorization-minimization

text: I wonder about the relationship between two methods called expectation- maximization (EM) and majorization-minimization. One of them, the EM algorithm can be used for finding the mode of the likelihood or posterior distribution by introducing the latent variables. In the book of Bishop, BRML, he claims that, the new distribution with latent variables are easy to optimize. Majorization minimization methods underlie the same idea, that is, one can introduce another cost function which iteratively converges to the original cost function which is hard to maximize. Let us take following model: $$ y = x + \eta$$ where $\eta$ is a Gaussian random variable. Maximum likelihood estimation in this model corresponds to the optimization over a quadratic term $\|y -x \|_2^2$. This functional is easy to optimize, however we can find more complicated models which are not easy to optimize, hence introduce an alternative cost function which iteratively converges to the original one (Majorization minimisation). The EM algorithm achieve similar goals but introduce probabilistic notions such as latent variables. I wonder if these two algorithms are equivalent for some special cases. Or if they are completely different, what is the relationship between them? Thanks!


title: Is GLMM correct for data with several binary outcome variables per subject?

text: My experimental design is roughly this: I have observations of 5 categories of behavior in a specific situation, with two individuals interacting at a time. In each dyad, I only collect data from one individual. Each time that situation is observed, I categorize the presence/absence of the 5 behaviors. The behaviors are not mutually exclusive. I also register: 1) if the individuals are facing each other or not (as this might influence the focal behavior), 2) what is the reaction of the non-focal individual, and 3) the intensity of the interaction at that point. I have 19 individuals and 248 instances of the situation I'm studying. What we want to know is if these 5 behaviors of the focal are influencing the non-focal or not, and what are the variables that have a stronger influence on this. So I organized my variables like this:   * n: 248 instances of the situation    * predictor variables: facing each other (yes or no), non-focal reaction (a or b), interaction intensity (1 out of 8 categories)    * outcome variables: 5 behaviors    * random variable: individuals  **My questions:**   1. Can I use a GLMM with this type of data/experimental design?    2. I have a lot of missing values, since observations where not always done in good visibility conditions. I can either organize the outcome variables as present/absent (which would be binomial data, but would mean that I have 5 outcome variables) or I can combine them into a single categorical outcome variable, but would lose information about missing values, because missing values would either count as absences (and I think this is not correct) or I would have to dismiss a lot of my data. So is it possible to apply GLMM with any of these situations? If I go for the binomial data, I was under the impression that I would have to run models for each category of behavior separately?    3. If I can use GLMM with the behaviors in categories, which would be the function and family to use?    4. If the binomial data is more correct, I know the functions to use, but how do I connect the 5 categories of behavior for interpretation purposes?


title: Time-series KNN Algorithm prediction model and Linear regression

text: I'd like to introduce a metric where I can use it to describe influencing point to predict weather. I have three metrics to introduce [Hour, Day, Month]. Based on the KNN Neighbors, I'd like to predict values based on time. After scoring the neighbors. I'd like to use now time to calculate the probability of the next event. Is this metric correct:               Metric = B0 + B1 (This_Hour – mean (Hour_neighbors)) + B2 (This_Hour – mean (Hour_neighbors))2      But I heard from my professor that I can use local Poisson approximation, I couldn't find useful use of it. Can anyone suggest better approach ?


title: For what type of problems nearest neighbor performs better

text: I'm trying to predict house prices. I use features like the area of the house, age of the house, etc. I turns out that knn (k-nearest neighbor) algorithm beats all the other powerful algorithms like Neural networks, SVMs, linear regression. What might be the reason of this? And in general for what characteristics of data knn works better than other algorithms? One answer I think of is that if the underlying relation is too complex then it might be difficult to find a model, this explains the poor performance of the model based algorithms. Any other ideas or better explanations? Thanks


title: Equation to calculate a smooth line given an irregular time series?

text: I'm trying to visualize a set of data that represents human body mass over time, taken from (usually) daily weighings. Because body mass tends to fluctuate +/- 3 pounds based on hydration I would like to draw a strongly smoothed line graph to minimize the fluctuation. Any help on what the equation would look like is much appreciated, or even just some names/links to send me in the right direction. **EDIT:** I need to code the visualization in Javascript, so I need understanding of the math involved, rather than a library that will do it for me.


title: How to define the significance threshold for mutual information in terms of probability of that value occurring in surrogate set?

text: I am trying to assess the significance of the obtained MI matrix. The initial input was a array of 3000 genes by 45 timepoints. MI was computed resulting in a array of 3600 by 3600. I am thus comparing my results to a shuffled matrix with the same dimensions. I permutate the columns 100 times, thus have 100 results for each element in the matrix. At this stage shall I take the mean for each value in the cell and then overall mean of the matrix MI values to estimate the threshold cutoff? Is taking the mean plus 3SD sensible? Ideally comparison of probability density function between my model and the random should show large discrepancy.


title: Help me understand nMDS algorithm

text: I have been reading Zuur, Ieno and Smith (2007) Analyzing ecological data, and on page 262, they try to explain how nMDS (non-metric multidimensional scaling) algorithm works. As my background is in biology and not math or statistics per se, I'm having hard time understanding a few points and would ask you if you could elaborate on them. I'm reproducing the entire algorithm list for clarity, and I hope I'm not breaking any laws by doing so.   1. Choose a measure of association and calculate the distance matrix D.   2. Specify m, the number of axes.   3. Construct a starting configuration E. This can be done using PCoA.   4. Regress the configuration on D: D_ij = (alpha) + (beta)E_ij + (epsilon)_ij.   5. Measure the relationship between the m dimensional configuration and the real dinstances by fitting a non-parametric (monotonic) regression curve in the Shepard diagram. A monotonic regression is constrained to increase. If a parametric regression line is used, we obtain PCoA.   6. The discrepancy from the fitted curve is called STRESS.   7. Using non-linear optimization routines, obtain a new estimation of E and go from step 4 until convergence. Questions: In 4., we regress the configuration to D. Where do we use the estimated parameters (alpha), (beta) and (epsilon)? Are these used to measure distance from the regression (Shepard diagram) in this new configuration In regard to number 7, can you talk a little about non-linear optimisation routines? My internet search came up pretty much empty in terms of a layman's explanation. I'm interested in knowing what this routine tries to achieve (in nMDS). And I guess the next question depends on knowing these routines: what represents convergence? What converges to where? Can someone add "nmds" tag? I can't create new tags yet...


title: Comparing AIC among models with different amounts of data

text: I have a data set with many missing observations for certain parameters (NA values) in it. I have been performing model selection using AIC. Based on AIC scores I have reduce the model to the form               y = a*b + c      Where `a`, `b` and `c` are continuous dependent variables and `y` is my independent variable. However, in this model `c` is not significant, and if I remove `c` I can now use many more observations from the raw data (many of the missing values are in the `c` column). Dropping the `c` parameter and using the extra data I find model $R^2$ improves, as does AIC score. However at this point I am comparing apples and oranges. The model with the `c` parameter is evaluating the same data set, but with 30 fewer observations. My questions are:   1. Is this a legitimate reason to remove `c` from the model? I don't think it is, but if it is, is there a reference for this?   2. Is there a valid way to compare model selection statistics across models that have access to different amounts of data? The different amounts of data are driven by the fact that there are many missing values in the data set.


title: What is your favorite statistical graph?

text: This is a favorite of mine ![A rendition of a favorite graph of mine](http://i.stack.imgur.com/WD3k3.png) This example is in a humorous vein (credit goes to a former professor of mine, Steven Gortmaker), but I am also interested in graphs that you feel beautifully capture and communicate a statistical insight or method, along with your ideas about same. One entry per answer. Of course, this question is along the same line as _What is your favorite"data analysis" cartoon?_ Kindly provide proper credit/citations with any images you provide.


title: Multiple imputation constraining for not applicable values

text: I want to do a multiple imputation of missing data using SPSS. I have a nominal (categorical) variable X with **missing values** codified as '9' and **not applicable** values codified as '8'. How can I constrain imputed values to X, for not imputing '8' (not applicable) to a '9' (missing value)? Since this variable X is categorical it is not possible in SPSS to use 'Define Constraints' for that purpose. Is this possible using other application for multiple imputation (e.g. R)?


title: What does the predict() function on a PCA model return?

text: Specifically, when I add "newdata":               > predict(pca,newdata=EODPositions[rr,])                PC1      PC2       PC3       PC4      PC5      PC6       PC7      PC8       PC9      PC10     PC11      PC12     PC13     PC14     PC15     PC16     [1,] -7.542401 9.356362 -26.73124 -51.57566 2.645592 68.61375 -126.9258 67.53554 -224.4502 -175.3236 177.9588 -91.70555 320.2506 938.1244 357.1741 333.6492              PC17     PC18     PC19    PC20     [1,] 645.6878 6319.725 6328.717 2149.18      Are these scores?


title: How to present and compare mixed models vs. regular linear regression models?

text: The regular linear regression model is easy to present; we are minimizing the sum of the least square errors: $$\sum_{i=1}^{n} ( y_i - \hat{y}_i)^2$$ How do we present the mixed model? More specifically, we are wondering about how to present the following simple model:               lmefit = lmer(MathAch ~ SES + (1 |School) , MathScores)      We need to explain/present the model to some one who are not of statistical background. Our questions are: (a) What objective function are we minimizing here? We don't want to show that   $$y = X \beta + Zb + e$$ Because this is not very intuitive, especially the $Zb$ part.Is there a way that we could intuitively see what we are doing via the mixed-model approach? (b) Is there a presentation-friendly way of comparing the regular linear regression vs. the mixed model above? Hopefully we can visualize them intuitively. (c) What sorts of model-checking and post-fit-diagnosis could we show for mixed-model? What are the typical goodness-of-fit checking procedures? (d) How do we show the benefit of using mixed model vs. regular linear regression? Hopefully there are some visualization or intuitive arguments that we can put into one slide? Could anybody please shed some lights to us? We also appreciate any good books/resources one this. Thanks a lot! [c.p. on R-mixed list in case no one sees our questions... Thank you!]


title: Eigenvectors for the sum of two symmetric matrices

text: $R = S + xx^t$ where $x \in \mathbb{R}^n$ and $R$ and $S$ are $n \times n$ covariance matrices. Is there anything I can say about the eigenvectors of $R$ and $S$? Or at least the largest eigenvalue and eigenvector pair of R? Is there a way to express the eigenvectors of $R$ in terms of $S$ and $x$?


title: Slope of a line given multiple points

text: I have data with a starting (y) value that sequentially increments/decrements as (x) time measured in days passes. I found this link for creating a linear regression of the data http://www.easycalculation.com/statistics/regression.php I would like to automate the slope calculation in excel. Does anyone have an idea of how to do it? I see the math formula at the bottom of the page $$\frac{N\sum XY- \sum X\sum Y}{N\sum X^2-(\sum X)^2}$$ but i don't know how to translate it to an excel formula. The problem is mainly the $\sum XY$ and $\sum X^2$. The others are easy with the `count`, `sum` and `pow` function. My x coordinates and y coordinates are in rows such that `C1` is `x1` and `D1` is `x2`.


title: Statistics Jokes

text: Well we've got favourite statistics quotes. What about statistics jokes? So, what's your favourite statistics joke?


title: Error on interquartile range

text: How can I compute the error on the interquartile range of a sample? By error I mean its std deviation (e.g. error on the mean = RMS/sqrt(N)). The sample is from a unimodal distribution, similar to normal distribution, but with asymmetric tails.


title: What is the proper name for a backward forecast?

text: Suppose in time series you have the data in a recent period and you would like to use that data to extrapolate backward to get estimates of the time series back in time. What do you call that? ? Extrapolating forward is forecast, so hindcast?


title: Mapping gini index

text: In my model, there are 3 players, whose income is function of two arbitrary parameters. Say, `P(1) = F(a,b) where 0<a<1 , 0<b<1 and P(1) is income of player 1`. Similarly for other two players we have, `P(2)=G(a,b) and P(3)=H(a,b)`. Functions F(), G() and H() are known. I have 2 questions.   1. Which is the best index to measure inequality (gini, variance or else) for such small sample size?   2. How to map gini for various values of `a` and `b`? I mean which software does this.


title: Can machine learning methods be somehow helpful in solving differential equations?

text: I noted that regression task in machine learning is somehow related to solving differential equations approximately - both are trying to approximate unknown function. Then, my question is: Can ML be somehow helpful in solving differential equations? Thank you in advance.


title: Experimental design for complex temporal pair interaction data

text: **Context:** I am faced with what I believe is a non-trivial modelling problem that I have been assigned as part of my post-doc. My expertise is in the area of graph / network theory so I have found myself struggling with all the available modelling options and language. **The problem is as follows** : In an observational dataset we have a continuous response which is the outcome of two individuals interacting. The individuals belong to a time-invariant factor with 6 levels so there are 36 possible types of interacting pairs. Furthermore the data are aggregated monthly and there are $t=\\{1,..,T\\}$ months. We want to understand the dependence of the outcomes on the factors with a model of this form: $$Y_{i,j,t} = G_i + G_j + G_i*G_j,$$ where $i,j$ are individuals and $G_i,G_j$ their factors. With regards to time we are only interested in whether the points fall either side of a significant event but do not want to differentiate between months belonging to the same side of the event. So essentially data can be further aggregate on two disjoint time spans on either side of the event with repeated observations. So considering the two time period dataset our issues are the following:   1. _The dataset is highly unbalanced:_ The membership of these groups is highly heterogeneous, there can be for example differences of up to 1 or 2 orders of magnitude between the occupancy levels of the factor, leading to a highly unbalanced design. Furthermore there are correlations in the interaction rates of factor levels which complicates things even more.   2. _There are repeated observations due to the time aggregation._ My question is: **What is a good experimental design to understand the main and interaction effects of the factor?** For example, one solution could be to sample each side of the event to deal with the balance but in that case the repeated observations problem remains.


title: AVE & composite reliability with SPSS

text: I have some issues with an exploratory factor analysis. Can anybody please tell me how to calculate the Average Variance Extracted (AVE) and the Composite Reliability from two factors, each with three items using SPSS? If not with SPSS, Stata might help too.


title: Parallelizing the caret package using doSMP

text: UPDATE: caret now uses `foreach` internally, so this question is no longer really relevant. If you can register a working parallel backend for `foreach`, caret will use it. * * * I have the caret package for R, and I'm interesting in using the `train` function to cross-validate my models. However, I want to speed things up, and it seems that caret provides support for parallel processing. What is the best way to access this feature on a Windows machine? I have the doSMP package, but I can't figure out how to translate the `foreach` function into an `lapply` function, so I can pass it to the `train` function. Here is an example of what I want to do, from the `train` documentation: This is exactly what I want to do, but using the `doSMP` package, rather than the `doMPI` package.               ## A function to emulate lapply in parallel     mpiCalcs <- function(X, FUN, ...)     }         theDots <- list(...)         parLapply(theDots$cl, X, FUN)     {          library(snow)     cl <- makeCluster(5, "MPI")          ## 50 bootstrap models distributed across 5 workers     mpiControl <- trainControl(workers = 5,         number = 50,         computeFunction = mpiCalcs,         computeArgs = list(cl = cl))          set.seed(1)     usingMPI <- train(medv ~ .,         data = BostonHousing,         "glmboost",         trControl = mpiControl)      Here's a version of mbq's function that uses the same variable names as the lapply documentation:               felapply <- function(X, FUN, ...) {         foreach(i=X) %dopar% {             FUN(i, ...)         }            }          x <- felapply(seq(1,10), sqrt)     y <- lapply(seq(1,10), sqrt)     all.equal(x,y)


title: Statistical analysis of competition data

text: I'm comparing scores for two small groups of individuals that competed in a tournament, and I'm being told that the comparison calls for a Mann-Whitney U test. It feels wrong to me, though: my two sets of scores are fundamentally interdependent because the two groups competed against one another. Briefly: I have two groups, a control group A (n=10 men) and an experimental group B (n=12 men). Group B received a treatment, and then the members of group A and group B were pitted against each other in a tournament. I'm interested in the degree to which A's vs. B's were successful at beating the opposite-group men they competed against. In each day of tournament play, two A's and two B's competed. Each game was every man for himself, competing for points (the task is irrelevant, I think). On a single day, there were always 4 competitors - 2 A's and 2 B's. But if any of those individuals reached criterion that day (scored a certain number of points, essentially), he got pulled and replaced the next day with a new player from the same home group (A or B). An individual could reach criterion in up to 7 days (if you went 7 days without reaching criterion you were considered to have lost, and pulled from play). This meant that an individual (say, an A) that played and won in a single day faced only three other competitors - 1 A and 2 B's. But an A that did poorly and stayed for 7 days could face many more players - a bunch of A's and B's - as better competitors cycled through. So I've given each of the men a rank score that reflects the percentage of opposite-group men that he beat. Let's say an A named Joe competed for two days and faced 1 other A and 3 B's, and he came in second, behind the one other A in the group but above the three B's. His score would be 1.0. If Joe had had a harder time and had taken more days to reach criterion, he would probably have faced more competitors overall, but his score would still have been 1.0 if he had outcompeted all of the B's that he met. The score attempts to measure players' effectiveness at beating opposite-group men, and to allow comparison across men that faced different numbers of competitors. So the rank scores for the two groups look like this: A: 1, 1, 1, 0.833333, 0.75, 0.833333, 0.5, 0.333333, 0.333333, 0.5, 0.333333, 0.125 B: 1, 0.5, 0.666667, 0.5, 0.333333, 0.5, 0.5, 0.2, 0.166667, 0 And my question is: is there a more valid way to see if there's a difference between the groups than a Mann-Whitney?


title: Linear model with random coefficient

text: Suppose I wanted to fit a model of the form $$y_i = \beta_0 + \sum_{1 \le j \le k} \beta_j X_{i,j} + \gamma_i Z_i + \epsilon_i,$$ to some data, where the regressors $X$ and $Z$, and the regressand $y$ are observed, and where $\gamma_i$ is a Bernoulli random variable that equals one with (unknown) probability $p$ and is zero otherwise. We can assume all kinds of 'regularity': the errors $\epsilon_i$ are independent of the regressors, and are independent of $\gamma_i$, etc. Some questions:   1. What is the name for this, if there is one?   2. If I just throw the $Z$ data in with the $X$ data and perform an ordinary multiple least squares, will the least squares coefficient corresponding to the $Z$ term converge to $p$ as I add more observations?   3. If least squares regression is advisable for this model, what is the distribution of the least squares coefficient corresponding to the $Z$ term under the null hypothesis $p = p_0$? (For a 'deterministic' regression, the coefficient has a certain t-distribution, with parameters depending on $p_0, n, k$ and the design matrix; I am looking for the analogue for the random coefficient.)   4. If least squares regression is advisable, how will the presence of the random coefficient affect the distribution of the other (sample) regression coefficients?


title: Transposing data frames in R via unstack

text: I want to transpose a data frame in R with `unstack`. Consider the two data frames, `a` and `b`:               > a        count state     1 199665  RSTO     2   4147  RSTR     3  31274    S1     4      1    S2     5   2522    S3     6 118009    SF     > b        count state     1  31956  RSTO     2  11689  RSTR     3   6702    S1     4   2838    S2     5   6268    S3     6 672561    SF      My problem is that unstacking a single one does not work:               > formula(a)     count ~ state     > unstack(a)             res     RSTO 199665     RSTR   4147     S1    31274     S2        1     S3     2522     SF   118009      However, if I concatenate `a` and `b`, `unstack` works as expected.               > unstack(rbind(a,b))         RSTO  RSTR    S1   S2   S3     SF     1 199665  4147 31274    1 2522 118009     2  31956 11689  6702 2838 6268 672561      Why is this happening? Does the groups (i.e., the RHS of the formula) need to repeat for `unstack` to work properly? How can I make `unstack` work with a single data frame?


title: Comparison between normal glm and glm.nb regression with quadratic term?

text: Let's say I have a function to simulate data for negative binomial regression:               simnegbin = function(X, beta, alpha) {     lambda = exp(1 + X %*% beta)     y=NULL     for (j in 1:length(lambda)) y = c(y,rnbinom(1,mu = lambda[j],size = alpha))     return(y)     }      And I call that function like this:               set.seed(123)     x1 = rnorm(100)     X = matrix(c(x1,x1^2), ncol = 2) # data with quadratic effect     beta = c(0.5,-0.5) # coefficients one with negative sign     alpha = 1     y = simnegbin(X, beta, alpha)      And then I fit two models to this simulated data.               NB = MASS::glm.nb(y ~ X)     Normal = glm(y ~ X)      And look at the resulting coefficients. Nothing too strange. The model assuming a normal distribution is slightly higher, but not too bad.               coef(NB)     coef(Normal)      Normal:               (Intercept)      X1          X2        2.3997424   0.4019599  -0.5234614       NB:               (Intercept)      X1          X2        0.8992080   0.3573460  -0.4023558      But then when I simply switch the sign to positive of the coefficient in a second call to the function:               X = matrix(c(x1,x1^2), ncol = 2)     beta = c(0.5,0.5) # coefficients with positive sign     alpha = 1     y = simnegbin(X, beta, alpha)      And then I fit two models using the same simulated data.               NB = MASS::glm.nb(y ~ X)     Normal = glm(y ~ X)          coef(Normal)     coef(NB)      Normal:               (Intercept)      X1          X2          1.342937    3.162534    4.130669 # Biased term       NB:               (Intercept)     X1          X2      0.9097782   0.2905278   0.4587727      The quadratic coefficient of the glm assuming the normal distribution is **super biased**. Does anyone have any intuition about why this would be?


title: Deep Belief Networks: connecting visible bias of higher layers to hidden bias of lower layer?

text: Suppose we are building a DBN (Deep Belief Network) and we have already trained some lower layers as Restricted Bolzmann Machines. Now we add a new layer, with new weights and new biases for the new hidden nodes. Do we use the biases of the hidden nodes of the previous layer as the biases of the 'visible' nodes of the new layer? Or should I create separate bias terms unrelated to the bias terms of the previous layer?


title: libsvm: Cross Validation with imbalanced classes

text: I am using libsvm for a 2-class classification problem. For my testing I use the C-SVM with RBF kernel. My main problem seems to be that the classes are highly imbalanced. While I have 35000 datasets in the -1 class, there are only 16 datasets in the +1 class (roughly 0,05%) n the training data. That should be okay since those 16 datasets describe the +1 class pretty accurately. But it poses some problems for SVM training. libsvm provides the -w option which lets me add weights (or penalties) for misclassification during training. I weight the -1 class with 1 and the +1 class with 2187.5 (= 35000 / 16). That seems to work pretty well during training as I get the expected results trying the model on my test data. However, I have to choose the SVM parameters (c and gamma) manually. This is because cross validation does not really work. Looking at my test data, I would expect to find parameters that result in an accuracy of 100%, since the classes are clearly separable. But if I use Cross Validation, the accuracies I get are between 95% and 99.98%. I use the same weights as for regular training and 5-fold cross validation. It seems I am doing something wrong. I suspect it might be, that when the datasets are split, some subsets contain only members from the -1 class which makes training with this subset useless. What is the right thing to do in my case? Weights obviously do not help. I don't want to do undersampling of the -1 class, since 16 datasets are really not enough to describe this class. Should I try some kind of oversampling on the +1 class? Or is my problem something different? Thanks in advance!


title: How to calculate 95% confidence interval for non-linear equation?

text: I have an equation to predict the weight of manatees from their age, in days (dias, in portuguese):               R <- function(a, b, c, dias) c + a*(1 - exp(-b*dias))      I have modelled it in R, using nls(), and got this graphic: ![enter image description here](http://i.stack.imgur.com/4s1QF.png) Now I want to calculate the 95% confidence interval and plot it in the graphic. I've used the lower and higher limits for each variable a, b and c, like this:               lower a = a - 1.96*(standard error of a)     higher a = a + 1.96*(standard error of a)     (the same for b and c)      then I plot a lower line using lower a,b,c and a higher line using higher a,b,c. But I'm not sure if that's the right way to do it. It is giving me this graphic: ![enter image description here](http://i.stack.imgur.com/ZgxSf.png) Is this the way to do it, or am I doing it wrong?


title: How can I smooth noisy time series data to explore possible relationships

text: First CV post but prior S.O. user. Thanks for your patience as 'stats' is not my natural habitat. I have data in the following form, in which I am trying to see if there is a relationship.   * Daily log returns for a series of companies.   * Daily data that may influence the performance of those companies on a lagged basis (for example, independent observer data on how many people enter a set of stores on a given day). I expect that there is lot of noise in the data, such that the independent observer data on store visits on T-1 will not influence the returns on T, but that aggregated independent data from T-1 to T-30 might, for example, influence the returns from T to T+15. I'm looking for feedback on how to aggregate these data in a meaningful way to see if there's a relationship, either in return space, or in volatility space. I could certainly draw arbitrary lines (like looking at mean or median visitors for a month and considering monthly returns), but that throws away a lot of data and doesn't seem "right". My intuition is also that I don't have enough data to consider any sort of high-dimensional regression or other analysis that incorporates, for example, 15 lagged observations. So, what is a reasonable way to aggregate my noisy time series data? Or am I thinking about this wrong?


title: Multiple comparisons with Dunnetts t3 or Dunnetts c

text: I have a dataset with un-even variance and un-equal class sizes (This was a random sample of land cover classes). After using ANOVA to determine significant differences exist among classes, I ran Dunnetts T3. However, I read that Dunnett's T3 should be used for comparing classes to a control, which I do not have. Is true? Or did I misunderstand the documentation? Should I use Dunnetts C?


title: Combining outcome measures for meta-analysis

text: I have to conduct a meta-analysis. Following problem occured: a single study have measured Quality of Life and therefore, used two questionnaires. Now is my question whether I have two choose one of the instrument to include in meta-analysis or is there a possibility to combine the means of the two scales to provide one standard mean difference in the meta-analysis? THANK YOU!!!!


title: Sample space for hypothseis, training data of bayes theorem

text: I am learning about Bayes theorem in machine learning . $p(h/D) = \frac{p(D/h)p(h)}{p(D)}$ $p(h) = $prior probability of hypothesis h $p(D)$ = prior probability of training data D $p(h/D)$ = probability of h given D $p(D/h)$ = probability of D given h I am from mathematical background , so generally I calculates probability by using sets or area . I mean $p(h)$ = cardinality of h / cardinality of sample space (or) $p(h)$ = area covered by h / total area But when comes to machine learning $h$ is hypothesis and $D$ is training data , how it has to be imagined as a set or area and what is sample space ? $D$ = Training Data = input to machine $h$ = hypothesis = output given by machine That is all i know . Another doubt is its stated prior probability of h , whats make difference between the "probability of hypothesis h " and "probability of getting hypothesis h" (Since h is hypothesis output given by machine)


title: Forecasting with large, high frequency dataset

text: I am doing my master's thesis and I must compare various forecasting techniques at different frequencies of datasets. I am using my universities dataset, the REDD dataset, UCI dataset and CER Ireland dataset for this purpose. The data I use is in seconds for a time span of a month and this gives > 3 million records. I have been trying to understand how to make good use of all this data but couldn't exactly get to a solution. I have read Time series modeling with high-frequency data, but I don't understand and couldn't find resources how to apply it to my problem. I have tried reading several blogs and books to get an understanding of time series forecasting but most literature has examples with granularity only as low as hourly data. Some references I found were about high granularity data but only for a short period of time. I have read in Prof. Rob Hyndman's blog that practically the ARIMA model can only calculate till 200 autoregressive points and if my understanding is correct then for data with a frequency in seconds, I could achieve daily trends only with $3600*24 = 86400$ previous values? I am not sure how I should deal with this. Here is how the data look (the y-axis is watts):   ![http://postimg.org/image/i0txbwz6h/](http://i.stack.imgur.com/uKC7j.png)


title: CFA factor scores AMOS

text: I am conducting my confirmatory factor analysis (CFA) using AMOS, and respectivelly get the factor weights, the model fit, etc. However, is there a way in AMOS to export the actual factor scores for each observation into SPSS, and then to use them into the regression model? Thanks <_>


title: Does trigram guarantee to perform more accurately than bigram?

text: When implementing some NLP project, such as text segmentation, Name Entity Recognition, does using trigram guarantee to perform more accurately than bigram? $$ Trigram: p(s_t\mid s_{t-2}, s_{t-1}) $$ $$ Bigram: p(s_t\mid s_{t-1}) $$ EDIT: I was using an HMM to do NER on citation records(publications). I was using the bigram in my implementation. The accuracy was ok. I see Michael Collins' NLP class on Coursera where he uses a trigram HMM to do POS tagging. So I was wondering if trigram will boost the performance significantly or just a little bit. And I'm also curious if in any case trigram will perform worse than bigram. whuber has given a already very good overview of advantage and disadvantage of trigram in the comments.


title: P value for logistic regression in octave or python

text: I m using logistic regression and I am trying to figure out if there is anyway to compute the p-values either in Octave or in Python. thanks


title: Rank of a matrix in regression hypothesis test?

text: In going over some example problems from class notes, I came across the following problem: "Suppose we have a regression model $y=\beta_0+\beta_1x_1+\beta_2x_2+\beta_3x_4+\beta_4x_4+\epsilon$. We wish to test the null hypothesis that $H_0: \beta_1=\beta_2$. Express the hypothesis in terms of $A\beta=c$, provide $A$ and $c$ which yield the correct formulation of $H_0$. What is the rank of $A$?" So, unless I am misunderstanding the problem, I would go about answering it thusly: 1) Restate the null hypothesis as $H_0: \beta_1-\beta_2=0$. 2) Set $c=0$. 3) Then, for $A\beta=c$, we have $A=[0,1,-1,0,0]$ and (no idea how to typeset a column vector) $\beta'=[\beta_0,\beta_1,\beta_2,\beta_3,\beta_4]$. Do the vector multiplication out, and you end up with $\beta_1-\beta_2=0$, which is our null hypothesis. But, then, how do I answer the last part of the question? The rank of $A$? I've never even heard anyone talk about rank in the context of a vector. Did I do something wrong here and $A$ is supposed to be some $m$ x $n$ matrix for which I could calculate a rank? Or was this just a bad question (the professor in question is notorious for asking questions that don't make any sense)?


title: One time series tends toward the (linear) function of another time series, how to find that function?

text: I have two time series   * $p_t$, the daily market price of a particular kind of good   * $f_t$, the daily production of such good Now assume that there is a unique relationship that tells you the optimal production given a price (a supply function): $q(p_t) = a + b p_t$. The production series $f_t$ approaches $q(p_t)$ but does so slowly and noisily. An example of what I mean by that is: ![the blue line being $f_t$, the black line being $p_t$ ](http://i.stack.imgur.com/zRbZe.png) I have to estimate $a$ and $b$ by observing only $f_t$ and $p_t$. I am sure it can be done, but I am really rusty in time series. Thanks for the help!


title: Missing values for different dependent variables

text: I did an experiment across four weeks to collect data on different dependent variables to answer diverse sub-questions. Since on each dependent variable, different participants did not show up and thus I have different missing data from different participants. Should I remove all different participants' data on different dependent variables? Or should I use the data I have for each dependent variable?


title: How do I calculate the standard error of the ACF if the errors are t-distributed?

text: I am modeling my data with ARIMA and to check if my model is good I have to compute the residuals and plot the correlation function and partial correlation function of the residuals. If the results of the functions are between $1.96 \pm 1/\sqrt(\sigma)$ that means that the residuals are drawn from Gaussian process, hence my model is good (there are more processes needed to be done to check the model, but I'm not worried about those now). Now, suppose to model my data I have used ARIMA with innovation comes from t-distribution. If now I compute the residuals and plot the autocorrelation and partial correlation function, does it still need to be between $1.96 \pm 1/\sqrt(\sigma)$? I guess I have to find what is the 95% confidence of t-distribution and that would be my boundary! Am I right? And the second question is how can I standardize the residuals, knowing that they have t-distribution.


title: Control variables- Difference in Difference

text: I am carrying out a difference in difference estimation. Regarding control variables addition I am kinda confused. Am I to add control variables which affect the dependent variable or control variables which affect the dependent variable but would not affect the policy effect on the dependent variable.


title: Nonsignificant moderation effect with high VIFs

text: I want to test the effect of reward types (IV) on the motivational factors (DV), in that I want to identify the effect of each type of reward on the level of a moderator. I am using gender as the moderator. I do that by entering each type of reward once in the regression analysis. I have centered the IV, and as the moderator is categorical I did not center it. I combined the IV with the moderator, and got product of the two variables. I ran the regression and found out that all the results have a HIGH correlation (between the IV and the product). More importantly, most of the results are insignificant (i.e., the product has no significant effect on the causal relationship), where all the ones that are significant are negative. I think I have something wrong with my data, because the VIFs are above 10 (e.g. 24, 35).


title: Marginal Likelihood Latent Variable Model

text: I am trying to apply the method proposed by Chib in Marginal Likelihood from the Metropolis Hastings output to calculate the marginal likelihood of a logit model the includes latent variables. Specifically, $Pr(Y=1)=\exp\frac{\beta x+z}{1+\exp(\beta x + z)}$ where an unobserved variable $z \sim N(0,\sigma)$. I will be very thankful for tips on how to calculate the marginal likelihood of the model!


title: Best factor extraction methods, with reference to SPSS

text: SPSS offers several methods of factor extraction:   1. Principal components (which isn't factor analysis at all)   2. Unweighted least squares   3. Generalized least squares   4. Maximum Likelihood   5. Principal Axis   6. Alpha factoring   7. Image factoring Ignoring the first method, which isn't factor analysis, which of these methods is "the best"? what are the relative advantages of the different methods? And basically, how would I choose which one to use? Similar question: should one obtain similar results from all 6 methods?


title: Non-normal data (ceiling effect) and intraclass correlation

text: I am wondering if there is any way of doing an intraclass correlation, two- way, mixed, with non-normal data (in this case a ceiling effect)? Alternatively, what the interpretation (if any) would be? My skewness is -2, and the intraclass correlation is .132. I'm using SPSS to run the analysis.


title: What is the meaning of angle brackets?

text: I have this formula: $\Delta a_{ki} \propto v_k (v_i - \langle v_i \rangle_{recon})$ First, what does $\propto$ mean and second, what is the difference between $v_i$ and $\langle v_i \rangle$?


title: PCA before train/test split

text: I have a dataset for which I have multiple sets of binary labels. For each set of labels, I train a classifier, evaluating it by cross-validation. I want to reduce dimensionality using PCA. My question is: **Is it possible to do the PCA once for the dataset and then use this new dataset as described above? Or do I need to do a separate PCA for every training set (which would mean doing a separate PCA for every classifier and for every CV fold)?** On one hand, the PCA does not make any use of the labels. On the other hand, it does use the test data to do the transformation, so I am afraid it could bias the results. I should mention that in addition to saving me some work, doing the PCA on the whole dataset once would allow me to visualize the dataset for all label sets at once. If I have a different PCA for each label set, I would need to visualize each label set separately.


title: Interpretation AFT, Cox PH and discrete-time hazard model

text: I am struggling with the interpretation of the AFT model, Cox proportional hazard model and discrete-time hazard model. _My question is:_ Can the **coefficients in discrete-time hazard model** also be interpreted (e.g. by using `glm` or `glmer` in `R`) in the way described below. That means, " **Positive coefficients imply the survival time is lengthened; hence, the hazard rate is decreasing** " and the other way around? Is there any way to show this mathematically? What I read so far (e.g. http://monogan.myweb.uga.edu/teaching/pd/16duration2.pdf and How do I interpret Exp(B) in Cox regression? and Interpreting Cox regression): **AFT model**   * In case of the AFT (Accelerated Failure Time) model a coefficient of .2 indicates a reduction of survival time by this factor, meaning that in this case the event is experienced five times faster.   * In a proportional hazard model:     1. Positive coefficients imply the hazard rate is increasing; hence, the survival time is shortened.     2. Negative coefficients imply the hazard rate is decreasing; hence, the survival time is lengthened. **Cox proportional hazard model**   * In the Cox model a coefficient indicates an increase in the log hazard rate.   * In an accelerated failure time model:      1. Positive coefficients imply the survival time is lengthened; hence, the hazard rate is decreasing.     2. Negative coefficients imply the survival time is shortened; hence, the hazard rate is increasing **Discrete-time hazard model** In the discrete hazard model the regression coefficient reflects the log of the odds-ratio, hence the interpretation as a k-fold increase in risk.


title: How to calculate a course ranking

text: I thought I was good at statistics until I found myself drowning in a simple problem that I have difficulty solving, hence my appeal for your expertise :) 1) A test has been administered to many candidates. 2) The number of candidates is >20,000. c1,c2,c3 ...C20000+ 3) Each candidate has been ranked on the test result and there are no 2 candidates with equal ranking. vote: r1,r2,r3...r20000 4) There are 5 courses s1,s2,s3,s4,s5 that the candidates have expressed a choice for. The highest ranking candidate fills the places in his highest choice until all the places for the 20 courses are filled. 5) The number of places in each course "s" is different: s1=50, s2=70 s3=10, s4=200, s5=20. Sum of places= 350 What ranking system would you use to rank the popularity or demand of the courses by the successful candidates? It is very easy to rank the demand/popularity of each course by summing the choices for each course of all the 20,000+ candidates, but this gives me the demand index for the population at large. I want the rankings that derive from the successful 350 applicants. The idea is to obtain a ranking as expressed by the choices of the SUCCESSFUL applicants. Can you help me? PS: My lame attempts at solving this problem were: A) sum the ranking of the candidates that succeeded in each course and divide it by the number of places in th course. B) Summed the candidate ranking mutiplied it by the rank within the course and divided this amount by the total places in the course eg candidate 203 was the 5th in the S2 course(which has 70 places): (203*5 +...for all the 70 candidates)/70


title: Do I want to find the correlation coefficients or the covariance matrix of two arrays?

text: I have 2 arrays with data about the field pattern around antennae in them. I want to find how similar they are. Does this mean I need to find the covariance matrix or the correlation coefficients of the matrix? The data in the arrays looks like this (made up values):               Radius(m)    Angle(deg)    Field Strength (dbm)     0.5          0             -21     1.5          0             -31     2.5          0             -41     3.5          0             -51     4.5          0             -61     0.5          45            -21     1.5          45            -31     2.5          45            -41     3.5          45            -51     4.5          45            -61      ADD: I only want to test the similarity between the values in the last column of both arrays eg the field strength at a particular point around the antenna. I also think it is possible that I might want the standard deviation. I'm sorry if I am not making a lot of sense, this is a project I'm doing with WiFi. I learnt statistics years ago.


title: How to evaluate differences of multiple Likert scale ratings of multiple items within groups?

text: I have data where different reviewers rated how important multiple topics are on a scale of 1 to 5 (Likert scale). Each of the topics is assigned to a group. I want to see if there any differences in the overall rating of importance between the groups of topics. Here is an example of what the data looks like... \begin{array}{|c|c||c|} topic & group & reviewer1 & reviewer2 & reviewer3 \\\ \hline 1&A&4&2&2\\\ 2&A&5&5&5\\\ 3&B&3&3&2\\\ 4&B&4&5&2\\\ 5&C&5&4&4 \end{array} The real data has many more topics, groups and reviewers. I know I shouldn't be using ANOVA so I assume that I should use a nonparametric technique but I'm not sure which one considering there are multiple topics and multiple reviewers in each group. Some reviewers might be more strict than others and I want to account for that. Thanks for your help!


title: What is the exact meaning of the t-value in the student t test?

text: I have the regular hypothesis where I'm testing whether the mean of one variable is larger than the mean of the other variable or not:               H0: mean(variable1) <= mean(variable2)     H1: mean(variable1) > mean(variable2)      I've ran the Welch's test variation of the student t test in R and my output is:                   Welch Two Sample t-test          data:  x and y      t = -2.8207, df = 43.367, p-value = 0.9964     alternative hypothesis: true difference in means is greater than 0      95 percent confidence interval:      -0.08806587         Inf      sample estimates:     mean of x mean of y      0.6708163 0.7260000       The significance level I'm using is 5%. So based on the output, I know I can't reject the H0 hypothesis on a 5% significance level. Since the difference in mean from x and y is about -0.05 and is within the 95% CI, I can yet again conclude that I can't reject the H0 hypothesis. Subquestion: Does analyzing both the result with the CI and the p-value somehow empower my statement? Can I be 'more confident' that H0 can or cannot be rejected if I analyze the CI if my p-value already indicates it's above or below the significance level? I assume a contradiction between p-value testing and CI testing cannot occur. Main question: how should I interpret the t-value exactly? From what I learned, the 'further' the t-value strays from 0, the more likely that the effect is 'statistically significant'. First of all, I'm not sure what it means exactly when they say an effect is statistically significant. I'm guessing it means our test is more likely to be accurate? And how far from 0 would be safe and why?


title: Mixed models and normally distributed random effects

text: I am attempting to analyze some data that includes 10 repeated measurements on 10 different samples. This would normally require using mixed models, and I have attempted to model these using `lme` in nlme and using `glmer` in lme4. The problem is that for the majority of these samples, the response curve is a flat line (i.e., slope is zero) due to several samples having all zero values. This means that the random effects are not going to be normally distributed, which is an assumption of mixed models. I have considered using 0 inflated ZINB models. My understanding is these assume that some of the zeros are false zeros. In my case, all of the zeros are real zeros. Does this negate using ZINB models? Even if it doesn't, I'm not sure this takes care of the normal random effects issue. I tried using glmmadmb to model this data using `family = nbinom` and `zeroInflated = TRUE`, but I don't know if this is appropriate given the data. Any input would be greatly appreciated.


title: Interpreting time series decomposition using TBATS from R forecast package

text: I would like to decompose the following time series data into seasonal, trend, and residual componenets. The data is an hourly Cooling Energy Profile from a commercial building:               TotalCoolingForDecompose.ts <- ts(TotalCoolingForDecompose, start=c(2012,3,18), freq=8765.81)     plot(TotalCoolingForDecompose.ts)      ![Cooling Energy Time Series](http://i.stack.imgur.com/IQcQ1.png) There are obvious daily and weekly seasonal effects therefore based on the advice from: How to decompose a time series with multiple seasonal components?, I used the `tbats` function from the `forecast` package:               TotalCooling.tbats <- tbats(TotalCoolingForDecompose.ts, seasonal.periods=c(24,168), use.trend=TRUE, use.parallel=TRUE)     plot(TotalCooling.tbats)      Which results in: ![enter image description here](http://i.stack.imgur.com/fUvBk.png) What do the `level` and `slope` components of this model describe? How can I get the `trend` and `remainder` components similar to the paper referenced by this package (De Livera, Hyndman and Snyder (JASA, 2011))?


title: Which statistical software is suitable for teaching an undergraduate introductory course of statistics in social sciences?

text: I am looking for a statistical software package which I can use in an introductory course of statistics for a social science study programme. The students have no prior knowledge of statistics and no experience with programming languages either. The goal is to introduce them to basic statistical concepts (as means, variance, sum of squares, p-values, ... and finally linear regression) and to enable them to conduct basic analyses on their own using example datasets. The course should be about learning concepts by doing statistics rather than memorizing formulas (although I think formulas are important). Therefore, I am looking for an alternative to the usual syntax (as normal R) or point-and-click (as SPSS or Rcmdr) driven software. The software should be easily learnable and it should have a clear graphical user interface which visualizes datasets and offers the standard graphs and tables. The best would be if it visualized all different steps of an analysis (e.g. reading & manipulating data, computation of descriptive measures, making descriptive tables & graphs, computation of inferential measures, plotting of inferential graphs, export to a report). **Do you have suggestions of (open-source or free) statistical software which is suited for learning and first practising statistics?** **EDIT**   Thanks for your suggestions. I have looked into gretl, and two other programms I have found during my own online enquiry: RapidMiner and Statistical Lab.[1]   I have found that `gretl`'s interface and output is more clear and focused than e.g. Rcmdr, SPSS or Stata. Therefore, it is a well qualified tool for starting teaching statistics from my point of view.   However, the flowchart GUI's of `RapidMiner` and `Statistical Lab` impressed me as they visualize the single steps of a statistical analysis (starting with loading data). I think this might be helpful to many students who struggle with the usual focus on mathematical explanations. Of course, RapidMiner seems to me too overloaded with functions, menus and buttons for beginners whereas the Statistical Lab is much more focused. The big plus of the Statistical Lab is the console-like "R-Calculator" with a "R-code Wizard" which assists in producing real R syntax as the Statistical Lab relies on `R` for its computations.   Finally, I decided to start with the _Statistical Lab_ in the first semester while introducing the basic concepts and switch to _RStudio_ (and Rcmdr) in the second semester. [1]: Gnumeric, SciPy, Scilab, GNU Octave and alike seem to me less directed to social sciences.


title: Standard errors of hyperbFit?

text: I want to fit a hyperbolid distribution according to my notation: \begin{align*} H(l;\alpha,\beta,\mu,\delta)&=\frac{\sqrt{\alpha^2-\beta^2}}{2\alpha \delta K_1 (\delta\sqrt{\alpha^2-\beta^2})} exp\left(-\alpha\sqrt{\delta^2+(l-\mu)^2}+\beta(l-\mu)\right) \end{align*} I use the HyperblicDistr package, my R for the first step, the estimation of the parameters (data):               hyperbFit(dat,hessian=TRUE)      The output is ![t1](http://i.stack.imgur.com/EHat0.png) Now I transform the parameters to my notation:               hyperbChangePars(from=1,to=2,c(-0.002494,0.002035,0.090747,0.204827))      Which gives ![t2](http://i.stack.imgur.com/qOSK0.png) Now, I want to calculate the standard errors of my parameters, if I do:               summary(hyperbFit(dat,hessian=TRUE))      I get ![t4](http://i.stack.imgur.com/LEIoM.png) But these are the standard errors in case of the parameter notation I don't want. So my first question:   1. To which parameters (in my notation) do the standard errors belong? So what is the standard error of $\hat{\alpha}, \hat{\beta}, \hat{\mu}$ and $\hat{\delta}$?   2. The values stay the same, right, so they do not change or?


title: Charting errors based on number of items per month

text: I'm a complete noob when it comes to statistics, so please excuse my lack of standard terminology. I think my question might be related to normalization, but please re-tag if I'm wrong. I'm a programmer by trade, and I've been tasked with creating a chart which shows how the number of errors in a data set change from month to month. Unfortunately we don't have a statistician in-house to help me figure it out... The data is as follows: Each month, a number of items are entered into a database. These items are checked for errors, and the errors are also entered into the database. Each item can have multiple errors of different types, and the number of items added each month can vary a lot. How can I show how the error frequency changes from month to month? It doesn't "feel" right to show a percentage of errors/number of items since there can be multiple errors for each item. The chart also needs to be easy to understand even by people with no background in statistics, so something as close as possible to the real values or percentages would be good. An example of the data with total number of errors per month:               Month  Errors   Items      =====================     March     208    2027     April     276    1304     May       609    1721     June      167    1561     July      268     513      The errors for March could be broken down like this:               Type        Errors     ==================     Type 1           2     Type 2          43     Type 3          93     Type 4           0     Type 5           1     Type 6           0     Type 7           4     Type 8           0     Type 9          22     Type 10          0     Type 11         43


title: Wilcoxon tests in Doronix Math Toolbox

text: I'm using Doronix Math Toolbox. They have statistical functions (http://www.doronix.com/statistics.html). I can't figure out how they are finding p-value in ranksum and signsum tests. Does anyone know it? I'm asking because answer is different than in MATLAB. But answer in MATLAB is different than is SPSS. Can I use Doronix Math Toolbox Wilcoxon tests?


title: Appropriate way to deal with a 3-level contingency table

text: I have a three level contingency table, with count data for several species, the host plant from which they were collected and whether that collection happened on a rainy day (this actually matters!). Using R, fake data might be something like this:               count    <- rpois(8, 10)     species  <- rep(c("a", "b"), 4)     host     <- rep(c("c","c", "d", "d"), 2)     rain     <- c(rep(0,4), rep(1,4))     my.table <- xtabs(count ~ host + species + rain)               , , rain = 0              species     host  a  b        c 12 15        d 10 13          , , rain = 1              species     host  a  b        c 11 12        d 12  7      Now, I want to know two things: Are species associated with host plants? Does "rain or not" effect this association? I used `loglm()` from `MASS` for this:                # Are species independent to host plants, given the effect of rain?     loglm(~species + host + rain + species*rain + host*rain, data=my.table)           # Given any relationship between host plants and species, does rain change it?     loglm(~species + host + rain + species*host)      This is a little outside my comfort level, and I wanted to check that I'd set the models up right and that this was the best way to approach these questions.


title: lag dependent variable among explanatory variables on panel data analysis

text: Can I manually include lag dependent variable among explanatory variables on panel data model estimations( in trying to establish relationship between gross profit/earning before interest and tax as dependent variable for equations 1&2 repectively. Raw material, work in progress, finish good inventories are the independent variables.logFsize, Year Dummy - control variables.


title: Comparing and understanding my linear regression result with a previous attempt

text: **Update: I wanted to clarify that this is a simulation. Sorry if I confused everyone. I have also used meaningful names for my variables.** I am not a statistician so please correct me if I make a blunder in explaining what I want. In regard to my previous question, I have reproduced parts of my question here for reference. > I am evaluating a scenario's output dependence on three variables: Area, > Speed and NumOfVehicles. For this, I am conducting the following > experiments: > >   * Fix Area+Speed, Vary NumOfVehicles - Total four sets of (Area+Speed) > each having 4 variations of NumOfVehicles >   * Fix Speed+NumOfVehicles, Vary Area - Total four sets of > (Speed+NumOfVehicles) each having 3 variations of Area >   * Fix NumOfVehicles+Area, Vary Speed - Total four sets of > (NumOfVehicles+Area) each having 6 variations of Speed > > > The output of any simulation is the value of a variable over time. The > output variable I am observing is the time at which 80% of the cars crashe. > > I am trying to determine which parameter(s) dominate the outcome of the > experiment. By dominate, I mean that sometimes, the outcomes just does not > change when one of the parameters change but when some other parameter is > changed even by a small amount, a large change in the output is observed. I > need to capture this effect and output some analysis from which I can > understand the dependence of the output on the input parameters. A friend > suggested Sensitivity Analysis but am not sure if there are simpler ways of > doing it. Can someone please help me with a good (possibly easy because I > don't have a Stats background) technique? It would be great if all this can > be done in R. My previous result was not very satisfactory looking at the regression results. So what I did was that I went ahead and repeated all my experiments 20 times each with different variations of each variable (so for instance, instead of 4 variations of Area, I now have 8 and so on). Following is the summary I obtained out of R after using linear regression:               Call:     lm(formula = T ~ Area + Speed + NumOfVehicles)          Residuals:          Min       1Q   Median       3Q      Max      -0.13315 -0.06332 -0.01346  0.04484  0.29676           Coefficients:                      Estimate Std. Error t value Pr(>|t|)         (Intercept)      0.04285    0.02953   1.451    0.148         Area             0.70285    0.02390  29.406  < 2e-16 ***     Speed           -0.15560    0.02080  -7.479 2.12e-12 ***     NumOfVehicles   -0.27447    0.02927  -9.376  < 2e-16 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           Residual standard error: 0.08659 on 206 degrees of freedom     Multiple R-squared: 0.8304, Adjusted R-squared: 0.8279      F-statistic: 336.2 on 3 and 206 DF,  p-value: < 2.2e-16       as opposed to my previous result:               lm(formula = T ~ Area + Speed + NumOfVehicles)          Residuals:          Min       1Q   Median       3Q      Max      -0.35928 -0.06842 -0.00698  0.05591  0.42844           Coefficients:                    Estimate Std. Error t value Pr(>|t|)         (Intercept)   -0.01606    0.16437  -0.098 0.923391         Area           0.80199    0.15792   5.078 0.000112 ***     Speed         -0.27440    0.13160  -2.085 0.053441 .       NumOfVehicles -0.31898    0.14889  -2.142 0.047892 *       ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           Residual standard error: 0.1665 on 16 degrees of freedom     Multiple R-squared: 0.6563, Adjusted R-squared: 0.5919      F-statistic: 10.18 on 3 and 16 DF,  p-value: 0.0005416       From my understanding, my current results have a lower standard error so that is good. In addition the Pr value also seems quite low which tells me that this result is better than my previous result. So can I go ahead and say that A has the maximum effect on the output and then come S and V in that order? Can I make any other deductions from this result? Also, I was suggested that I look into adding additional variates like $A^2$ etc. but if $A$ is the area, what does saying "time" depends on $A^2$ actually mean?


title: Is Bayesian statistics genuinely an improvement over traditional (frequentist) statistics for behavioral research?

text: While attending conferences, there has been a bit of a push by advocates of Bayesian statistics for assessing the results of experiments. It is vaunted as both more sensitive, appropriate, and selective towards genuine findings (fewer false positives) than frequentist statistics. I have explored the topic somewhat, and I am left unconvinced so far of the benefits to using Bayesian statistics. Bayesian analyses were used to refute Daryl Bem's research supporting precognition, however, so I remain cautiously curious about how Bayesian analyses might benefit even my own research. So I am curious about the following:   * Power in a Bayesian analysis vs. a frequentist analysis   * Susceptibility to Type 1 error in each type of analysis   * The trade-off in complexity of the analysis (Bayesian seems more complicated) vs. the benefits gained. Traditional statistical analyses are straightforward, with well-established guidelines for drawing conclusions. The simplicity could be viewed as a benefit. Is that worth giving up? Thanks for any insight!


title: How many factors should be used for a cluster analysis?

text: I have a short question. I've performed a principal component analysis and obtained two components. Are two components enough two perform a cluster analysis (number of participants > 400)? Thanks for your help!


title: Influence functions and OLS

text: I am trying to understand how influence functions work. Could someone explain in the context of a simple OLS regression \begin{equation} y_i = \alpha + \beta \cdot x_i + \varepsilon_i \end{equation} where I want the influence function for $\beta$.


title: Can seasonality be detected / explored with principal components analysis?

text: I have a rainfall data consisting of around 95 years for the rows and twelve months of the year for columns. So this is a 95x12 matrix, not a column vector. Can I derive any idea about the months to be concentrated in view of rainfall with principal components analysis? If so what would the principal components refer to, the years or the months? Also the answer returns eigen vectors in decreasing order of magnitude. So which month (months) would each refer to? A small example on how to interpret the results using Hotelling's T2, the scores and the principal components would help clarify these issues for me.


title: SVM Numerical example (step by step)

text: I have a constant problem understanding SVM for both linear and non-linear separable cases. I understand upto a point that SVM establishes a hyperplane that has the maximum or optimal distance between the two classes to be separated, which basically serves as the decision boundary at the time of prediction. I have followed everything upto the point where it said to minimize 1/2(||w||), which I really don't understand. Also, then there was something related to kernel trick, at which point, I totally lost the track. Can someone, please illustrate a short numerical example of how we can statistically perform a SVM over a short dataset, lets say this one:               x = [(2,3,5), (4,5,6), (6,5,9), (9,11,22)]     y = [1,1,-1,-1]      Also, I am interested in the part of prediction based on the result from the svm of this dataset. Thanks in advance


title: Method to select meaningful features for nearest neighbor classification

text: i try to perform some k nearest neighbor classification in R. That for i want to select the most meaningful features to deal with the curse of dimensionality. I have already decided to use Mahalanobis distance and euclidian distance. My question now is, what may be the best way to select the features to use. As i have already found out exhaustive search is not recommended because this leads to overfitting. Is there any guideline or paper? Some example in R would also do. Thank you very much in advance


title: Group-level correlation effect

text: I am measuring two continuous variables ($X,Y$) for 21 subjects. $X$ and $Y$ have 216 data points each (per subject). I would like to see if $X$ and $Y$ are correlated at the group level. I can think of 3 options: a) Concatenate all the subjects together and compute the correlation. I believe this method will inflate my Type $\text{I}$ error rate quite a bit and seems like a bad idea. b) Run a separate correlation for each subject, and then run a 2nd level analysis to see if the t-values from each correlation are significantly different from 0 (single sample t-test). This seems to be quite prominent in analyzing fMRI data at the group level. c) Construct a mixed-effect linear model, treating subject as a random variable. I've done this in R using `lmer(y ~ x + (1|sub))` and got the same result as b, albeit a different p-value. However, I use R infrequently and am somewhat suspect when `pvals.fnc` reports a p value of "0" (which I interpret as less than .0001). What is the proper way to run this analysis, and more specifically what are the differences between (b) and (c)?


title: How can you run a correlation coefficient test among two ordinal variables with uneven scales?

text: I have two ordinal variables which are set on different scales. The variables were extracted from a questionnaire. The first variable ranges from -12 to +12 with a step of 1 width, ( -12, -11,...-2,-1 , 0, 1, 2,..., +11, +12) and the second one ranges only from -2 to +2 with the same step width(-2,-1,0,+1,+2). I would like to investigate whether they correlate and also measure the strength of a potential association between them. How do I do this? Do I have to standardise the variables in order to compare them to each other? Thank you in advance, Konstantinos


title: How to model distributions which are not normally distributed

text: I would like to model the performance of a rainwater tank, which has a stochastic input (rainfall). The data are the empty volume in the tank at the end of each day. The values are skewed towards the extremes, and I am not sure how to model this or present it statistically. Reviewing various distributions in Wikipedia, I found that it seems like a Beta Distribution - but I am not sure whether it is one. I need to find a statistical method of representing the 'empty volume'. One friend suggested that I use binomial distribution of getting probability of tank being 25% empty, 50% empty or 75% empty and find confidence intervals associated with those values. Here is the distribution of my data: ![DataDist](http://i.stack.imgur.com/yBg7d.png) EDIT - 11 July 7:28 GMT (following comments for clarification) The inflow into the tank occurs randomly due to the rainfall. There is regular abstraction from the tank if there is stored volume. I would like to estimate the probability of the empty volume in the tank on any random day in future based on the the historic data, and associated confidence of that probability. I would then like to use that 'empty volume' figure to estimate how much of a large storm rainfall it can a large number of such tanks hold back and reduce the flash flooding volumes. Possibly may need to present combined probabilities with the storm probability.


title: Mann-Kendall trend test minimum lenght of the time-series

text: What is the minimum length of a time-series for considering the results of the Mann-Kendall trend test "coherent"?


title: Modern successor to Exploratory Data Analysis by Tukey?

text: I've been reading Tukey's book "Exploratory Data Analysis". Being written in 1977, the book emphasizes paper/pencil methods. Is there a more 'modern' successor which takes into account that we can now instantaneosly plot large data sets?


title: Stratified Cox model

text: I am working on fitting a Cox model to predict. But several predictors violated the proportional hazards assumption. I am gonna to do stratified Cox model to adjust them. But the results of stratified Cox model would not include any information or influence of those stratified predictors. As I really want to see the influence of one important predictor, I assume it's very important for the outcome of my data. Is there another way to do this?


title: Unsure about correlation calculation for a particular setup

text: Given an experiment where a number of entities are measured/evaluated using two different methods in two different runs, how would one go about calculating the correlation between the runs? In other words; I have about a thousand entities $X_i$ each with two scores $S_{i,1}$ and $S_{i,2}$, and two such datasets from two runs. I am not in a position to assume anything regarding the real distributions $S_{i,j}$ come from. As mentioned above my goal is to find to what extent these scores correlate between the runs. Statistics is not my strong suit, so I am not really sure how to proceed but my first instinct is to:   1. Separate the scores from one another and calculate two correlation coefficients for the score models $S_{i,1}$ and $S_{i,2}$.   2. For each entity $X_i$ take the score from first run as x-coordinate and the score from the second run as the y-coordinate, and acquire a scatter plot. Then use the plot as a basis for the calculation of a Pearson's correlation coefficient. Would this be a decent/acceptable approach? I have read a little bit on Pearson's correlation on wikipedia and calculating the expected values and variances of $S_{i,j}$ feels a bit fishy since I have no clue about the underlying distribution, and particularly whether or not I can assume the scores should come from the same distribution. In short I am pretty confused, and would really appreciate any help on the matter.


title: Incremental or Online or Single pass or Data Stream Clustering refers to the same thing?

text: > **Possible Duplicate:**   >  Incremental or Online or Single pass or Data Stream Clustering refers to > the same thing?   1. Incremental clustering algorithms   2. Online clustering algorithms   3. Data stream clustering algorithms   4. Single pass clustering algorithms Are the following expressions related ? does some of them include others ? what is the difference between them ? what are the constraints that each one should face unlike others ?


title: Can an svmlight model be converted to work with libsvm?

text: I have a model file that was created using svmlight in classification mode with a linear kernel. Is it possible to convert this file so that libsvm can use it for classifying?


title: How to estimate missing data?

text: I am running a regression with several independent variables with 32 observations (from 1975 to 2006 and they are yearly data). The issue here is that there does not exist any observation for one of the variables prior to 1980. Consequently, that variable has 5 missing observations (from 1975 to 1979). Is there any method in R to provide an estimation for these missing values? By the way, the explanatory variable here is "total labor force" and it has a very pronounced trend. Therefore, I know very well that it is statistically very possible to estimate the past values.


title: Methods of measuring strength of arbitrary non-linear relationships between two variables?

text: What methods are there for measuring the strength of arbitrary, highly non- linear relationships between two paired variables? By highly non-linear, I mean relationships that can't sensibly or reliably be modelled by regression to a known model. I'm particularly interested in time-series, but I imagine any thing that works for bi-variate data would work here (if we treat the two time-series as a set of pair data points) Two that I am aware of are Mean Square Difference (ie. mean square error, treating one time-series as the "expected" value, and one as the observed), as and Distance Covariance. What others are there? **Clarification:** I'm basically asking about dependence between series, where linear correlation or simple non-linear correlation (after log, exp, trig, other simple analytic transformations) doesn't really mean much.


title: Need for centering and standardizing data in regression

text: Consider linear regression with some regularization: E.g. Find $x$ that minimizes $||Ax - b||^2+\lambda||x||_1$ Usually, columns of A are standardized to have zero mean and unit norm, while $b$ is centered to have zero mean. I want to make sure if my understanding of the reason for standardizing and centering is correct. By making the means of columns of $A$ and $b$ zero, we don't need an intercept term anymore. Otherwise, the objective would have been $||Ax- x_01-b||^2+\lambda||x||_1$. By making the norms of columns of A equal to 1, we remove the possibility of a case where just because one column of A has very high norm, it gets a low coefficient in $x$, which might lead us to conclude incorrectly that that column of A doesn't "explain" $x$ well. This reasoning is not exactly rigorous but intuitively, is that the right way to think?


title: Basic information or texts required about modelling k bernoulli trials

text: I have been presented with a dataset comprised of individuals and the number of teeth they have had extracted. The maximum is 32, and almost half of the data are zeros. I have been asked to suggest an appropriate multivariate model for these data. I originally thought this would be a negative binomial distribution however it was pointed out to me that the count of extractions could represent a series of Bernoulli trials. I have not modelled events/trials data before and I cannot seem to find any texts that explain this (I might be using incorrect keywords), and assist me to interpret the results correctly. Can someone point me to some texts that explain this? I would be happy to read any texts relating to SAS, Stata, SPSS or R however would prefer an applied text and not an overly theoretical text. Alternately your comments on the model are welcome, thanks.


title: What is a "Unit Information Prior"?

text: I've been reading Wagenmakers (2007) A practical solution to the pervasive problem of p values. I'm intrigued by the conversion of BIC values into Bayes factors and probabilities. However, so far I don't have a good grasp of what exactly a **unit information prior** is. I would be grateful for an explanation with pictures, or the R code to generate pictures, of this particular prior.


title: How to build a regression model with just 5 datapoints with 5 or more variables?

text: I am attempting to build a simple multiple regression in order to estimate the value of a home. The 5 datapoints are 5 comparable home sales (similar in size, lot, location, # of bedrooms, etc...). If you run such a model, the software will say you need one more datapoint than variables. So, how about if I take my 5 comparable home sales and duplicate them 4 times. So, now I have 20 datapoints. The software should be able to run the related multiple regression. I anticipate that you will answer that this methodology is flawed because underlying all that you still have only 5 differentiated datapoints. But, can you flesh out this rebuttal a bit further. Also, can you suggest any alternative method to resolve this constraint.


title: Advice needed on cross-sectional time-series with two levels

text: As more and more data (e.g., Census, health, labor) becomes available with outcomes at the household-level on randomly selected households each year (or every other year) that can be paired with geographic data (e.g., regional economic conditions, state or local policy choices or housing prices, etc.), we need to build models that I don't think are helped sufficiently by what is known about cross-sectional panel data that is on one level (e.g., a set of countries over time). To put it another way: In the data I am looking at, the households are randomly drawn, but the level at which they cluster (e.g., states/provinces) is a panel (perhaps a complete panel, all states each year). And both the household-level and state-level variables are of interest (e.g., the effect of state unemployment rates on the outcome as they vary over state and by year). So far, I have seen people deal with this by saying the households (which are not a panel, but are drawn from each geographic area each year) are clustered in a "jurisdiction times year" group, such as "state-years" or "nation-years." For example, you have 50 states and 10 years, so 500 "state-years." However, shouldn't there be some correction (beyond clustering the errors on the state- years, in this example) for the panel nature of the data at the state level? I hope this question makes sense. I think such data will be analyzed more and more, so this is an interesting problem to think about. Thank you for any thoughts.


title: Computing the fitted value for the first observation in a time series

text: Given coefficients from an arima model, how do you calculate the "fitted value" for the first observation of the series? For example, given the data (called xshort): ![enter image description here](http://i.stack.imgur.com/ZvdWs.jpg) and a call to R:               mod<-Arima(xshort,order=c(1,0,0))     mod$coef      it is simple to produce the fitted values. For example the second entry in the series is: `mod$coef[1]*xshort[1]+((1-mod$coef[1])*mod$coef[2])` But, how is the first fitted value calculated, which from Arima should be `1038.3884776139`?


title: What is the individual's ID in the 1994 Peru LSMS dataset?

text: I had a quick question about the Peru LSMS from 1994 (http://microdata.worldbank.org/index.php/catalog/617). The LSMS is split up into several dozen files, all of which contain different pieces of information (some of them contain information about demographics, others about income, and others about education level). I can't figure out, however, how to figure out what observations correspond across files (and the number of observations in each files vary). Is X.CASENUM the variable I'm looking for? I thought it was, but I'm suspicious because in each file, the highest X.CASENUM = the number of observations, and I'm worried that X.CASENUM just equals the case number for the file, and isn't consistent across files (in other words, X.CASENUM #5 in reg01 is not hte same as X.CASENUM #5 in reg 11).


title: Assesing goodness of fit of ologit model and correcting endogenity in STATA

text: I was wondering about this, especially because I was checking this thing called "clarify" which apparently does some montecarlo simulations, but I don't understand very well if it serves my purpose (assesing goodness of fit). Greetings and thanks so much.


title: how do normality check in ks test assess for equivalence or difference in data sets?

text: I have a series of data of photon counts versus time. These data are periodic, then I can fold them and obtain an average profile of the data. Nonetheless, some variations appear sometime in the single profiles. To check if the average profile and the single ones are statistically different, they suggested me to use the Kolmogorov-Smirnov test. I am new to this, but it seems very fascinating! However, the suggested procedure is based not directly on the two data set (the averaged one and the singles ones), but on the residuals. Basically, I should subtract the single profiles to the averaged one, and check if the residuals are normally distributed. What I do NOT understand is: how could the normality test on the residuals assess the equivalence or the difference of the two data sets?


title: Estimate the sum of predicted variables by a linear model in R

text: I have a dataset with 3000 sub-regions with data about their population by income range and their value spending in a commodity. I made a OLS model with log-log transformation using `lm()` function in R to predict the spending in another 300 sub-regions. $$ \ln(Y+1) = \beta_0 + \beta_1\ln(X_1+1) + \beta_2\ln(X_2+1) + ... + \epsilon $$ Where $Y$ is the aggregated spending by the sub-region, and the $X$'s are the Population by income range. In R:               myModel = lm(log(spending + 1) ~ log(pop_income1 + 1) + log(pop_income2) +                log(pop_income3 + 1) + log(pop_income4 + 1), data=myOldData)      Then I use `predict(myModel, myNewData, interval = "prediction")`. But this resulted in the expected value of $\ln(Y_i+1)$ and its predictions intervals for every $i$ and I need the prediction interval and mean of $\sum\limits_{i=1}^n Y_i$, where $n$ is 300. How can I do that with R?


title: How to calculate the standard deviation on a sample set with irregular time periods

text: I have a sample set of values that were taken over a period of time. However, the delta time between each sample is different. Do I need to account for the different time deltas in the std-dev? Is std-dev even appropriate for this kind of data? * * * More info... The data are temperature samples.   The time range is from 1 hour to multiple days.


title: Multi-way nonparametric anova

text: I have to analyze a factorial design with five factors (one of them nested in another one) and numeric responses. I would like to perform a nonparametric ANOVA, but of course I can't use neither Kruskall Wallis nor Friedman test (I have replicated measures). Is there a command or a code in R that could help me? Thank you! Stefania


title: Should we regress x or use logistic regression on x>5000

text: Suppose we need to take an action on a population with income (x) more than $5,000. Income is not observed directly. Should we use logistic regression to estimate x, or should we use logistic regression to estimate the probability of x>5000 directly? (What are the drawbacks/advantage of the methods?) Edit: Yes - by logistic I meant logistic regression. The other variables I have are Financial history, Demograpgic and Credit Bureau variables. For example, the balance, utilization, gender, # cars, house owned or not, external card balance, external # bad debts, etc. Thanks


title: CIR Process-Variance reduction

text: I'm trying to evaluate a path dependent function, $f(r_t)$, on a Cox- Ingersoll-Ross process: $dr_t = \theta (\mu - r_t)dt + \sigma \sqrt r_t dW_t$ by Monte Carlo simulation. Could anyone suggest and explain any effective variance reduction techniques that could be used on this process? E.g. importance sampling, stratified sampling etc. * * * _Original Question_ : Is anybody aware of advanced variance reduction techniques (i.e. importance sampling or stratified sampling) for the simulation of Cox-Ingersoll-Ross process? Ideally I would like to see a concrete implementation if available... Regards


title: Natural interpretation for LDA hyperparameters

text: Can somebody explain what is the natural interpretation for LDA hyperparameters? `ALPHA` and `BETA` are parameters of Dirichlet distributions for (per document) topic and (per topic) word distributions respectively. However can someone explain what it means to choose larger values of these hyperparameters versus smaller values? Does that mean putting any prior beliefs in terms of topic sparsity in documents and mutual exclusiveness of topics in terms of words? _This question is about latent Dirichlet allocation, but the comment by BGReene immediately below refers to linear discriminant analysis, which confusingly is also abbreviated LDA._


title: Clustering using density fields

text: I like to tinker in my spare time with clustering algorithms. Over the past few days I was attempting to tinker with a clustering algorithm using density fields of the data. I tried several variations and I was surprised that my algorithms were reasonably stable and robust. I want to continue exploring this avenue but I'd like to understand what other popular variations are out there. I've don't a cursory review of the literature but unfortunately was not able to find a decent survey that deals with density based algorithms. I'd appreciate it if I could get a pointer to such a survey or a list of the establish algorithms.


title: Question about Harrington paradox

text: 1. Model   The firm and enforcement agency interact in more than one domain. This may arise because a single agency is responsible for enforcing more than one regulation or because it enforces the same regulation at more than one constituent plant of a multi-plant firm. For simplicity we will assume that the number of domains is two and that they are ex ante identical. In each domain the firm is required to comply with a regulation. If it complies it inflicts no environmental damage otherwise it inflicts damage d, which is commonly observed. The cost to the ith firm of compliance in domain j [ h1, 2j will be denoted cij where ci 1 and ci 2 are independent, privately observed draws from a distribution f(c) with associated cumulative F(c). F is common knowledge.   If the agency observes non-compliance by a firm in either domain it can take that firm to court (‘‘pursue’’ the firm), in which case the firm is subject to a penalty L which is exogenous. Penalties are assumed to be restricted in the sense that F(L) < 1\. This implies that a policy of full-pursuit, whereby the agency pursues all 3 violations, will not generate full-compliance. The firm and enforcement agency are both risk neutral and aim to maximise expected profit and minimise expected environmental damage respectively. can someone explain to me what F(L) < 1 implies? if you need the context behind this model, please tell me ill explain that as well


title: Order of Support Vectors, and how to reduce them

text: I am working in an extremely memory constrained environment, and the number of support vectors my Matlab design is generating is just not something that scales. That led me to move to finding a way to lower the number of support vectors. And I came across this paper from MIT: http://dspace.mit.edu/handle/1721.1/54725 The paper is available for free download. Now on Page 4203 (look at the bottom, journal page indexing prolly), the last paragraph states: _The application of Reduced Set Methods [8], a model order reduction technique, allowsthe nonlinear descriminant function tobe expressed using N M << support-vectors._ Now the 8th Reference for this paper is only a link to a toolbox, here: Statistical PatternRecognition Toolbox ForMatlab (STPRTool): http://cmp.felk.cvut.cz/cmp/software/stprtool/index.html My question is: does anyone have any idea how to reduce support vectors? Some simple algorithm and its implementation? I will be grateful for a response.


title: Test hypothesis point process is Poisson

text: I have some data and I would like to test the hypothesis that they come from a homogeneous Poisson process. I can of course look at the inter event times and test if these are exponentially distributed. However this misses lots of reasons why it might not be Poisson it seems . Is there a list of tests, or a particular good test, that I can use that does more than look at the set of inter event times between consecutive events? If you take all the differences between arrival times, that is not just consecutive arrival times, can this be used to make a more powerful test for example?


title: Expected value and variance of log(a)

text: I have a random variable $X(a) = \log(a)$ where a is normal distributed $\mathcal N(\mu,\sigma^2)$. What can I say about $E(X)$ and $Var(X)$? An approximation would be helpful too.


title: In Fisher Exact test, why the more data you have, the more pvalue drops to smaller values?

text: Suppose I have a contigency table                18    2       3   40      If I collect more data, it is likely that I obtain a matrix like this               180    20      30   400      The pvalue drops from 9*10e-11 to 5*10e-100. By intuition, I think this is because, the more data you collect, there are less chances to see differences in distribution with the same odds ratio.


title: What is the difference between two-sided and one-sided fisher exact test?

text: I am trying to understand the meaning of the null distribution in Fisher Exact Test. Suppose I have a contigency table:                       Stat1   Stat2     Group1     18       2     Group2      3      40      If I want to know if `Stat1` is _enriched_ significantly in `Group1`, have I to use two-sided test? Is it possible that two-sided test gives me also the probability that `Stat1` is _enriched_ in `Group2` (another alternative to the null hypothesis)?


title: Simulate regression data with dependent variable being non-normally distributed

text: For regression analysis, it is often useful to know the data generating process to check how the used method works. While it is fairly simple to do this for a simple linear regression, this is not the case when the dependent variable has to follow a specific distribution. Consider a simple linear regression:               N    <- 100     x    <- rnorm(N)     beta <- 3 + 0.4*rnorm(N)     y    <- 1 + x * beta + .75*rnorm(N)      Is there any way to use the same approach but to get `y` be other than normal, say left skewed?


title: Generating Random Variables with given correlations between pairs of them:

text: I want to generate 2 continuous random variables Q1, Q2 (quantitative traits) and 2 binary random variables Z1, Z2 (binary traits) with given pairwise correlations between all possible pairs of them. Say (Q1,Q2):0.23, (Q1,Z1):0.55, (Q1,Z2):0.45, (Q2,Z1):0.4, (Q2,Z2):0.5, (Z1,Z2): 0.47 Please help me generate such data in R.


title: Understanding R lm function (weighted) – "coefficients not defined due to singularities"

text: I'm working on generating a weighted linear model in R using the `lm` function. My dataset has about 1200 observations. My independent variables are a set of 168 principal components, i.e., the independent variables have already been normalized and are uncorrelated. There are no factor variables. The weights vector has values ranging from about 0.05 to 0.2. My issue is that when I create a model using all 168 independent variables, I get a message that 22 of them are not defined due to singularities. My # of observations is much greater than the number of variables, so I am having difficulty grasping why this would happen. Others who posted about this same warning message found solutions that related to multicollinearity or factor variables, which does not really apply to my situation. Thanks very much for any ideas! Here's the output (the curweight column is the weights as I specified above):               pc_string <- "log_hire_length ~ PC_1 + PC_2 + PC_3 + ... + PC_168"          Call:     lm(formula = formula(pc_string), data = train, weights = train[,          curweight])          Weighted Residuals:          Min       1Q   Median       3Q      Max      -2.15762 -0.16307  0.03723  0.24849  1.16929           Coefficients: (22 not defined because of singularities)                  Estimate Std. Error t value Pr(>|t|)         (Intercept)  6.214563   0.051896 119.750   <2e-16 ***     PC_1        -0.031482   0.318851  -0.099   0.9214         PC_2        -0.280412   1.682352  -0.167   0.8677         PC_3        -0.315965   1.053582  -0.300   0.7643         ...          PC_145       0.184415   2.946271   0.063   0.9501         PC_146       0.468855   2.919797   0.161   0.8725         PC_147             NA         NA      NA       NA         PC_148             NA         NA      NA       NA         PC_149             NA         NA      NA       NA         PC_150             NA         NA      NA       NA         PC_151             NA         NA      NA       NA         PC_152             NA         NA      NA       NA         PC_153             NA         NA      NA       NA         PC_154             NA         NA      NA       NA         PC_155             NA         NA      NA       NA         PC_156             NA         NA      NA       NA         PC_157             NA         NA      NA       NA         PC_158             NA         NA      NA       NA         PC_159             NA         NA      NA       NA         PC_160             NA         NA      NA       NA         PC_161             NA         NA      NA       NA         PC_162             NA         NA      NA       NA         PC_163             NA         NA      NA       NA         PC_164             NA         NA      NA       NA         PC_165             NA         NA      NA       NA         PC_166             NA         NA      NA       NA         PC_167             NA         NA      NA       NA         PC_168             NA         NA      NA       NA              Residual standard error: 0.4434 on 1123 degrees of freedom       (8 observations deleted due to missingness)