title: Does it matter how you pull marbles out of this vat?

text: There's this large, well-mixed vat containing some unknown but finite varieties of marbles: $$\{v_{1},v_{2},v_{3},...,v_{k}\}$$ Some varieties are more common than others. Marbles are going to be pulled out of this vat by some machine. Each time you run the machine, it produces a set containing $q$ different kinds of marbles: $$\{v_{q_{1}},v_{q_{2}},v_{q_3},...,v_{q}\}$$ We're going to look at a bunch of these sets with the idea of estimating how likely it is that any $v_{n}$ shows up in such a set. So now we have to go and collect our data. The thing is, we don't know precisely how the machine picks marbles; it could pick up $q$ simultaneously, then discard and pick again if they aren't all different. Alternatively, it could keep picking up marbles until there were $q$ varieties. It could also pick a marble, then ignore the rest of it's kind. There's plenty of ways it could be picking the marbles. (Assume the vat is large enough that any finite number of marbles we pull out doesn't affect the population to a measurable degree.) The only thing we know about how the machine picks is that it doesn't discriminate against the marbles it picks as long as they are relatively unique. In other words, the likelihood that a given variety shows up in the $q$ set is determined by it's commonness in the vat, not the machine. This brings me to my question, will the results of our statistical analysis be affected at all by the manner in which the machine collects marbles?


title: What is the most effective way to populate the PCA training matrix for eigenfaces?

text: I am experimenting with the eigenfaces algorithm, and I am unclear on a number of the finer points of the technique. For starters, consider the matrix of images used to do the initial PCA. Is it "better" to have only images of different individuals in that matrix? Or is it useful to include multiple images of the same individual? Some of the demo implementations I have seen appear to be using only one image per individual. Some use more than one (and generally these seem to have N sample images per individual rather than varying sample sizes). The paper describing the algorithm seems to suggest that multiple images could produce more accurate results. It would intuitively seem like multiple images of the same person would tend to weaken the PCA results. For example, if I had 100 images where 90 were of the same person, the principal components would be more about distinguishing between those 90 images of the same person than distinguishing her from the other 10 individuals. Right? Does that same reasoning still apply when there is a uniform number of sample images per individual? If not, is there a guideline for how to balance the number of samples per individual against the total number of training images?


title: How to compute standard deviation of difference between two data sets?

text: I am running an experiment where I collect two data sets and I wish to measure the difference between the two. The two data sets are independent, with unknown probability distribution, and may not always have the same length. Calculating the mean difference is easy as pie, but i also want a measure of the standard deviation and I'm not sure how to go about doing that. At first I just got the std deviation of the second data set minus the average of the first set, but in retrospect, I'm not sure that is entirely correct. Any advice out there? As an example of the data:               3.98   4.39   4.09   4.31   3.81   3.67   3.94   3.90   4.39   3.60   3.99   3.53   3.82      vs               3.95   4.51   4.49   4.43   4.55   4.41   4.68   4.22   4.45   4.59   4.42      **Edit** : I want to etch this in stone a bit more to explain why ttest would not be the answer I'm looking for: All computations are done in Matlab. note: `var(one) = 3.61e-6, var(two) = 5.01e-06`.   Using gui11aume's answer:                std = 2.93e-3      when doing                [h,p,ci,stats] =ttest2(one, two)          stats.sd = 2.08e-3      realizing the variances are likely unequal, one should rewrite this to:               [h,p,ci,stats] =ttest2(one,two,0.05,'both','unequal')         stats.sd =  [1.9e-3 2.23e-3]      once again, thank you for your time


title: Contingency tables, Chi-Squared analysis and next steps in R

text: I'm hoping one of the more experienced R programmers/statisticians out there can help me with something, or explain why what I'm trying to do is stupid :-) In the past, I've done analysis of survey data in SPSS. I'll often create contingency tables that uses a computed segmentation variable to see if members of the segments exhibit different behavior along an outcome variable. Often the outcome variable is a response along a 5 or 7 point scale. Sometimes it might be a binary yes or no. Here is what my analysis would typically look like in SPSS: ![enter image description here](http://i.stack.imgur.com/oXaKF.png) In this case, I grouped respondents into one of 3 segments. I then re-coded a 5-point response to q39_1 into a binary grouping of 4/5 and 1-3. I did this because my n's were small in a few of the boxes. You can see that the response does seem to be dependent on the respondent's segment. The Chi-square test confirms this. I've managed to recreate this in R using the survey package with the following code:               > data <- read.spss("sample_data.sav", to.data.frame=T)     > design <- svydesign(ids=~1, weights=~Weight, data=data)     > svytable(~q39_1_recode2+alt_segment, design) -> table2     > prop.table(table2,2) -> table2prop     > table2                  alt_segment     q39_1_recode2         1         2         3                 1  68.05988  88.23581  40.40518                 2  31.97007 115.05503 168.27403     > table2prop                  alt_segment     q39_1_recode2         1         2         3                 1 0.6803950 0.4340373 0.1936234                 2 0.3196050 0.5659627 0.8063766     > chisq.test(round(table2,digits=0))              Pearson's Chi-squared test          data:  round(table2, digits = 0)     X-squared = 71.4848, df = 2, p-value = 3.001e-16      Now comes the snag. In the SPSS output, you'll notice that in the contingency table there are little subscripts that indicate where the column proportions are different from each other at an alpha of .05. I can toggle this column proportion on or off through a checkbox when I create the table. You can set the alpha level and then toggle a setting to adjust the p-values for multiple comparisons using the Bonferroni method. I'm obviously not a statistician, but my assumption was that this is essentially doing the z-test equivalent of pairwise t-tests that you might see from an ANOVA analysis. I find this type of analysis useful because not only do I now know that there is a relationship between the variables, I can now describe what that relationship is. **So, the question is: How do I replicate these pairwise z-tests in R while correcting for the multiple comparisons? Is this in fact the correct "next- step" after a chi-squared test of independence?** I've had trouble finding material that discusses what the next steps in an analysis should be after rejecting the null hypothesis in chi-sq test. I feel like I'm either missing something that's pretty easy or I'm trying to do something that I shouldn't be. Any help would be appreciated.


title: Deterministic algorithms inspired by MCMC?

text: There's a heuristic, which I find appealing, that says that for every stochastic algorithm, there should be at least one deterministic algorithm that performs better, provided the universe isn't adversarial. I first heard this principle articulated by Eliezer Yudkowsky here. Markov Chain Monte Carlo algorithms, which are some of the post powerful and general algorithms for approximating probability distributions are, seem like a counterexample to this principle. Has there been much investigation into algorithms which take their inspiration from MCMC algorithms (I'm mostly thinking of Metropolis Hastings and Hamiltonian Monte Carlo), but which are totally deterministic? If not, why not?


title: How to compare several learning algorithms in a single domain

text: I want to compare several learning algorithms in a single domain. In the book " Evaluating Learning Algorithms" I saw the procedures to compare two learning algorithms in a single domain or multiple algorithms in multiple domains. I assume that the procedure to evaluate multiple learning algorithms in a single domain would be the same as the procedure to compare two learning algorithms as long as I correct the error due to multiple comparisons. I would like to know how to perform this correction (for example Holm-Bonferroni) in non parametric test such as Wilcoxon or McNemar. For Wilcoxon, I think it that it would be no problem since I have the p-value, but how can I do it in the McNemar test? Any suggestions for an alternative method to compare multiple learning algorithms in a single domain would be appreciated as well. One possibility I can think of is ROC comparison as discussed in:   * DeLong E.R., DeLong D.M., & Clarke-Pearson D.L. (1988). Comparing the Areas Under Two or More Correlated Receiver Operating Characteristic Curves: A Nonparametric Approach. Biometrics, 44:837-45


title: Comparison of average values of data sets

text: I am working with two data sets of unequal sample sizes of 998 and 857. Average of first (998 samples), come out to be higher than the other dataset. To my surprise, when I split my complete data into unequal halves; first into 803 and 195 samples and other also into 819 and 38 samples. then comparison of average values of 803 subset of first complete dataset with the average of 819 subset of second complete data set showed a reverse trend in their mean values. Same reverse trend was observed with other subset of both the data set. My question is that is it possible that if mean of total items in A>B, their two subsets showed reverse trend in their means i.e. mean of A1 Is this because of unequal sample sizes???? or because of sample density distribution trend????? or both???? If this is possible also, then is there any way to explain this quantitatively?????? It would be really helpful, if anyone can help me on this..


title: Measuring representativeness of a sample using covariates

text: I was provided with quite a small sample of labeled (variable of interest) observations to train a model to predict unlabeled observations. All the observations are associated with many covariates. I'm assuming that the trained model will do better in measure of how well does this small labeled sample represent the unlabeled cases. Using only the covariates is there any way to measure if an unlabeled case will be poorly predicted or not? I can imagine that if the covariates are standardized and you measure the euclidean distance (covariates are continuous) from a unlabeled point to the labeled ones and this unlabeled one tends to be "far away" from the labeled sample the prediction accuracy would drop. I'm not really sure. Does anyone have any comments on this or any pointers on what to read to assess this or if my ideas are just totally off? By the way, the techniques I tried out for this task are Random Forests, MARS and boosted regression trees.


title: construct naive bayes classifier using defined characteristic

text: I found a paper about using naive bayes classification for identifying interested user. In that paper, they use some characteristic to define whether user is interested or not. However, they do not explain the step by step method to build the classifier. This is the link where i found the paper : http://ijcsi.org/papers/IJCSI-9-1-2-381-387.pdf Could anyone explain to me how they build a **naive bayes classifier** based on characteristic that they define? Do they build their own training data manually? I already try to contact the author but it seems still no answer.


title: Mixed Model t Test (Toothaker’s t-Test)

text: I'm trying to do a follow up test for a mixed model ANOVA, and I figured that I need do a Mixed Model t test also known as Toothaker's test. I just don't know how do it. I was wondering if anyone has an Excel sheet that can calculate this for me, because I know SPSS can't. I'm doing a mixed model ANOVA and I found significant interactions between my within subjects variables and my between subjects variables. Notes from a previous year's stats course told me to use a mixed model t-test to examine these interactions, and we were given an Excel sheet to help us calculate it. However I lost this Excel sheet. To clarify more I have two between subject conditions (drug/placebo) and two within subject conditions (memory tested at 2 different times). In my notes the test is called Toothaker's Mixed Model t-test and is described as such "Used when 1 (or more) within-subjects factor and 1 (or more) between subjects factor: Toothaker’s mixed-model t test (pools MS error)" Wikipedia seems to describe a similar thing "when there is a significant interaction between a between-subject factor and a within-subject factor, statisticians often recommended pooling the between-subject and within-subject MSerror terms.[2][page needed][citation needed] This can be calculated in the following way: MSWCELL = SSBSError + SSWSError / dfBSError + dfWSError"


title: Help interpreting Pearson's R p-values

text: I am trying to understand what these results mean for my analysis. The p-value seems extremely high, although the actual correlations do seem to reflect my data.               Pearson's product-moment correlation         data:  domain.of.choice_stats and objects.of.choice_stats      t = -0.4177, df = 10, p-value = 0.685     alternative hypothesis: true correlation is not equal to 0        95 percent confidence interval:     -0.6555880  0.4789388           sample estimates:               cor      -0.1309586      Here are the actual correlation results:                          Consumer.products    Services       Public.services            Gov        -0.95399809         -0.6622662      0.8703883                  Com         0.98559387          0.3404654     -0.9777812                Org         0.03533326          0.9271726      0.1740777              Here is the 'head' of domain.of.choice_stats and objects.of.choice_stats:                               Gov      Com      Org            Consumer        0        12        0                  Commercial      0        10        0                Health          0        7         3     Political       2        0         0                                 consumer.products    services   public.services            Consumer        10                    2           0                  Commercial      10                    2           0                Health          7                     6           1     Political       0                     0           2


title: Orthogonality of parameter space

text: I have a basic information theory question. I am fitting a highly parameterized model to some data. In general: $$ y = \sum\limits_{i=1}^{13} \alpha_i X_i $$ Currently I use gradient descent to find $\alpha$. The general assumption would be that the $X_i$s are i.i.d. However I am not sure if this necessarily holds true for my features and therefore need a method to analyse the correlation between $\alpha_i$s. Also if the feature space has high correlation what would be principled approaches to reduce model complexity?


title: Is there any test for a null hypothesis of non-normality?

text: I'm currently looking for a test having for null hypothesis that the sample does **not** come from observing a normally distributed random variable. In other words, I'd like to know if there's a test that's opposite to the K-S, Anderson-Darling, Shapiro-Wilk or Jarque-Bera.


title: SVM regression with LIBSVM and python - Expected runtime

text: Im new to regression, but quit experienced in classification and machine learning. In Classification, the state of the art technique for classification is SVM, so the first solution I would like to use for my regression problem is an SVM one. The regression problem im trying to solve can be described as follows: I have one dependent variable which i would like to estimate with 100,1000,4000 independent variables respectively, and estimate the accuracy with each size. I have around 100,000 observations which Im planning to use in a 90%-10% cross validation scheme.   The generating process is probably not linear, this is why I want to use SVM regression (I will test linear kernel, though) Question is, before I even attempt such process, what is the expected run-time for one regression iteration with linear or RBF kernel? Can anybody give a raw estimation? Is this problem feasible in the regression world?


title: How to do Univariate Heteroscedasticity Test

text: I just wanted to know how to do **_Heteroscedasticity Test on a Univariate Model_**?   * ex: an univariate autoregressive model   * ex: an univariate ARCH/GARCH model If it is possible, how does one do that in `R`?


title: What is the normal approximation of the multinomial distribution?

text: If there are multiple possible approximations, I'm looking for the most basic one.


title: Expected first hitting time for compound Poisson process?

text: Is there an analytic method to compute the expected first hitting time of an IID jump process with exponentially distributed time steps and jump sizes given by a Poisson random variable multiplied by a random variable that takes values {1,-1} with probability 0.5?


title: Comparing linear models

text: My Linear Models professor said that it is possible sometimes to obtain two different linear models for the same data, such that the models cannot be compared to determine which is better. How can this be? As an example he mentioned models that have different scales. What does it mean for models to have different scales, and why does it make them incomparable? Is scale mismatch the only reason for models to be incomparable? Is it only possible to compare a nested model to the full model that it derives from? If not, how does one compare two models that are not hierarchically related?


title: Age Interactions in mixed models

text: We are using a mixed model to describe trends as subjects age. We have multiple measurements per subject of a continuous outcome. We also adjust for the age at entry into the study. It is of interest if early entry and long follow-up are associated with different outcomes than late entry and short follow-up (i.e. same current age, but entered study at different ages) and so we are including an interaction. So our model is of the form: $y = current$ $age + baseline$ $age + current$ $age * baseline$ $age$ Are there any problems with estimation or prediction with a model of this form? I think that this is very similar to the standard age-period-cohort problem, but not sure if that means there are other considerations we need to have when modeling or interpretation.


title: What is the role of a categorical predictor in polynomial regression?

text: I understand that there is a function in R called `poly()` that can generate orthogonal polynomials--useful for applying on input variables before running a predictive model. My question is that what is the role of categorical variables when we generate polynomials? Are they to be excluded? ## Update: Dan, Thank you for your kind response. I'm not sure I understand it completely - let me explain the query in more detail. I'm trying to run logistic regression using glmnet on Titanic dataset.   Let us assume shortened set of columns: * class(factor with three levels 1, 2 ,3),    * sex(factor: male, female),    * Age (integer),    *survived(factor & target variable 0 or 1). The questions is it meaningful to create polynomial features based on these factors? e.g. class. If yes could you pls explain what it means?   I've seen examples with numeric input variables, where one can pass the entire input set to the poly() function and get polynomial features as output.   Your response is highly appreciated.


title: significant differences between groups in gender

text: I have a cell array containing the gender e.g.               GENDER = {'Male';'Female';'Male';'Male';'Male';'Female';'Female';'Female';'Male';'Male'};      and another array containing the respective groups e.g.               GROUPS = [1; 2; 1; 1; 3; 3; 1; 2; 2; 3];      I would like to find whether there are significant differences between the 3 groups in the number of male and number of female. Or maybe, even better, the differences between the ration male/female in the 4 groups. I think that the GENDER variable could be assumed Bernoulli-distributed so using the ANOVA should be statistically incorrect. If instead of the GENDER I've to analyse the difference in eye colour what is the best approach?               EYE_COLOR = {'Blue';'Blue';'Brown';'Black';'Black';'Blue';'Brown';'Black';'Black';'Black'};      Can I assume the EYE_COLOUR a categorical distribution? Thanks a lot


title: Training a Support Vector Machine in Matlab - number iterations exceeded

text: I'm trying to train a support vector machine on a set of training data (which is seperate to my testing data). `training_data` has size `MxN` where `M = 10175` and `N = 102` Each observation has a label `1` or `0`. (118 observations have label 1, 10057 observations have label 0). Initially, I ran the Matlab Statistics Toolbox method `svmtrain` like so:               svmtrain(training_data,labels,'Kernel_Function','rbf');      These are the default parameters recommended by Matlab. However after a while I get this message               ---------------------------     Error using seqminopt>seqminoptImpl (line 198)     No convergence achieved within maximum number of iterations.          Error in seqminopt (line 81)     [alphas offset] = seqminoptImpl(data, targetLabels, ...          Error in svmtrain (line 499)     [alpha, bias] = seqminopt(training, groupIndex, ...      _What next steps should I take?_ I tried using the `'Kernel_Function','linear'` option but the same problem happened. I'm not sure whether If I should use my own values for `rbf_sigma` or even for `boxconstraint`? If anyone can help guide me in getting a SVM to be trained then that would be really helpful! PS. Could there be something about my data that means the svmtrain function will _always_ exceed the maximum number of iterations, i.e. never converge? I don't want to aimlessly try different combos `svmtrain` method parameters, there are simply too many!


title: Ordered Response Variable

text: For regression with ordered response variable, there are different methods, for example, discriminant analysis, probit or logit model. I am wondering what are the different focuses of the different methods and which one is more often used.


title: Finding parameters to maximize expected utility of random variable

text: I'm trying to analyze some data consisting of five randomized parameters and a utility function which indirectly depends on the parameters, by experimentation. That is to say, the parameters of the experiment are chosen randomly, and successes and failures are counted up. I want to find parameters for which the expected utility of successes and failures is highest. From my days in calculus, I can see that an algorithm could consist of:   1. Regression to a (hopefully analytically tractable) surface   2. Finding a maximum   3. Finding the pre-image of my maximum (if I use any of the C libraries I've seen, which seem to focus on the maximum value, not its pre-image) But I'm not sure about the "fiddly bits" like:   * The distribution of points (I don't have any data yet)   * Any substantive idea of the shape of the surface, though I am expecting diminishing marginal utility, so it should be non-linear and have a bump.   * Numerical stability This seems like it should be straight-forward, in terms of applied decision theory. So, is my plan sensible? Any pointers to literature, algorithms, C or Haskell libraries? Addition in response to comment: I'm trying to find the "best" parameters in terms of student performance. The 5-tuple represents:   1. $b$: The "base" waiting time before seeing a problem again.   2. $p_1$: A constant factor if the student says the problem was easy.   3. $p_2$: A constant factor if the student says it was hard.   4. $p_3$: A constant factor if the student says it was "normal".   5. $p_4$: A constant factor if the student got it wrong. The waiting time for the next viewing is computed by multiplying all of the responses the student has issued, and the base waiting time, and computing $e^{b \prod p_{i,j}}$. So, for example, a wrong answer makes the waiting time much shorter. An 'easy' report makes it quite a bit longer. Now, if the student gets the next viewing wrong, we want to count it as a failure. If the student gets it right (regardless of the difficulty the student reports), we count it as a success. I want to maximize the utility function $\frac{|\text{success}|}{|\text{total}|}$ by varying the 5-tuple. I guess $\frac{|\text{success}|}{|\text{failure}|}$ would serve the same purpose.


title: goodness of fit in (prospective) relative risk regression

text: Instead of logistic regression when there is a prospectively follow cohort experiencing events there are methods to estimate a relative risk = RR (rather than an odds ratio) eg http://aje.oxfordjournals.org/content/159/7/702.full.pdf+html usually binary regression with a log link or a poisson regression this is imprtant as OR != RR when the event isn't rare and being prospective we know the prevalence so can estimate a proper RR. I wondered 1) whether there are any goodness of fit/residuals for these methods and 2) what the advantages over the COx model are - is it that these RR regression are good when the onset of the event isn't sudden so the actual event time isn't clear (though maybe interval censoring could be used) many thanks


title: Combining results with different confidence levels

text: I have an experiment where I throw a biased coin were with confidence interval [0.38, 0.42] and confidence level of 0.95 I will get heads. If I get heads then I throw a biased 6-sided die which gives me 1 with probability given in the interval [0.48, 0.52] and confidence level 0.90. What is the probability of getting 1 and what is my confidence level for that result. If no intervals and no confidence levels the problem is trivial and could be done with simple multiplication, but I am pretty sure simple multiplication will not be a solution in interval case.


title: difference between mixed effect and a simple ANCOVA

text: So I have collected data for 50 years of my dependent variable (continuous) data (x), independent (continuous) variable (y) across 10 locations. My null hypothesis is that the regression slopes between x and y is homogeneous across 10 locations (i.e. slopes of y against x is similar across all 10 locations) against the alternate hypothesis that slopes are heterogeneous (i.e. slopes of y against x is different for different States). I thought I could just do this by doing:                mdl1<-y ~ x*location       # includes interaction x and location      mdl2<-y ~ x + location     # no interaction between x and location      anova(mdl1,mdl2)      If the model with the interaction term is a better fit, then I can deduce that at least two locations have different slopes. But I really want to know which slopes are different for which I need to do a posthocs test. In this case what posthocs test should I do? Also someone told me that this analysis is wrong and that I have to use a mixed effect models with locations as a random factor. I tried reading about the mixed effect models but could not understand it. What I want to know is whether the analysis above is wrong and do I actually need to use a mixed effect model?


title: What are the correct edge case values of precision and recall and how to integrate them into cross validation?

text: I'm doing 10-fold cross validation on a dataset. But in some folds there are edge cases that the denominator in precision-recall calculation is zero (tp + fp =0). What are the correct values for precision and recall in this case? And what is the correct way of doing cross-validation (should I include these results when reporting average precision-recall over 10 folds)? PS: This question is very similar to What are correct values for precision and recall when the denominators equal 0? Thanks in advance : -)


title: What kind of statistical method is recommended to compare ratio in pre-test and post-test

text: I want to compare two different ratios. When this ratio is close to 1 it shows the better stability in the body. Would you please guide me that which statistical test is better for this situation? Please consider SPSS for answering this question. Here is some information about my situation: My study was conducted on 30 women. 15 women were takes as the experimental group and 15 women as the control group. Pre- and post-participation levels of core stability were measured by McGill core-stability tests. These tests were performed in four positions as flexory, extensory, right and left. Then a core-stability exercise program was performed for 8 weeks just for experimental group. After doing the exercise program four mentioned tests were taken again. Then the ratio of flexory/extensory and right/left were calculated. As much as this ratios are close to 1 it shows the better balance of core-stability. I know that comparing mean measure is not a right choice since for example 0.1 and 1.9 both shows bad balance of body but the mean of these values is 1 and it shows a perfect core-stability balance.


title: Binary logistic regression: low event rate

text: I want to investigate how people (lets say smokers and non-smokers) differ in various characteristics. So I want to perform a binary logistic regression with SPSS. However, my event rate is very low (number of events = 27 and number of non-event = 200). I want to include several predictors: 3 factors, which I gained from a PCA, 4 socio-demographic variables and 3 other variables (including categorial variables). I've read that if the sample size is too small, I can't perform a binary logistic regression.   * In my case, is it better to perform an exact logistic regression or a Firth's logistic regression due to the small sample size?    * If I would only consider the 3 components (without including the other 6 variables) in my binary logistic regression model , can I perform a binary logistic regression then? So Is the sample size here adequate? In that case, should I use the "enter" method or a stepwise method like "backward LR"? (In the theory model (backward LR), 2 out of 3 factors (p=0.030 (factor 1) and p=0.088 (factor 2) and p=0.000 (constant)) were included. In the overall model (enter), one factor was significant (p=0.33 (factor 1) and p= 0.94 (factor 2) and p=0.115 (factor 3) and p=0.000 (constant) . Which method should I prefer? Thanks for your help!


title: How to measure change in time series?

text: I have some weekly data about Apple app store ranking of the top apps listed (also of the sub categories such as games, weather, etc.). For example:                              week1   week2   week3   week4   ...     Angry Birds    15      13      1       5       ...     Weather Pro    ...      I am thinking about trying to estimate a parameter of "turbulence" or "change" for each category and compare them to each other. So, how can I measure change in time series?


title: using neighbor information in imputing data or find off-data (in R)

text: I have dataset with assumption that nearest neighbors are best predictors. Just a perfect example of two-way gradient visualized- ![enter image description here](http://i.stack.imgur.com/7frrm.jpg) Suppose we have case where few values are missing, we can easily predict based on neighbors and trend. ![enter image description here](http://i.stack.imgur.com/o2o1w.jpg) Corresponding data matrix in R (dummy example for workout):               miss.mat <- matrix (c(5:11, 6:10, NA,12, 7:13, 8:14, 9:12, NA, 14:15, 10:16),ncol=7, byrow = TRUE)     miss.mat          [,1] [,2] [,3] [,4] [,5] [,6] [,7]     [1,]    5    6    7    8    9   10   11     [2,]    6    7    8    9   10   NA   12     [3,]    7    8    9   10   11   12   13     [4,]    8    9   10   11   12   13   14     [5,]    9   10   11   12   NA   14   15     [6,]   10   11   12   13   14   15   16      Notes: (1) The property of **missing values** is assumed to be **random** , it can happen anywhere. (2) All data points are from single variable, but their value are assumed to be influenced by `neighbors` in row and column adjacent to them. So **_position in matrix is important_** and may be considered as other variable. My hope in some situations I can predict some off-values (may be mistakes) and correct bias (just example, lets generate such error in the dummy data) :               > mat2 <- matrix (c(4:10, 5, 16, 7, 11, 9:11, 6:12, 7:13, 8:14, 9:13, 4,15, 10:11, 2, 13:16),ncol=7, byrow = TRUE)     > mat2              [,1] [,2] [,3] [,4] [,5] [,6] [,7]     [1,]    4    5    6    7    8    9   10     [2,]    5   16    7   11    9   10   11     [3,]    6    7    8    9   10   11   12     [4,]    7    8    9   10   11   12   13     [5,]    8    9   10   11   12   13   14     [6,]    9   10   11   12   13    4   15     [7,]   10   11    2   13   14   15   16      The above examples are just illustration (may be answered visually) but the real example may be more confusing. I am looking if there is robust method to do such analysis. I think this should be possible. What would be suitable method to perform this type of analysis ? any R program / package suggestions to do this type of analysis ? ![enter image description here](http://i.stack.imgur.com/F0Mwv.jpg)


title: How to simulate with given probability?

text: I have a presumably simple question. What does this mean: with probability 0.3 variable X is simulated as > X~N(0,1). How can I simulate X using this?


title: randomness and computers

text: How do computers achieve randomness if they really do ? How to decide whether something is random or not? Is there a measure of randomness?                                     Thanks in advance..


title: How to add a non-negative constraint to lasso4j?

text: Lasso4j is a Java implementation of the Lasso L1-constrained fitting for linear regression. I would like to add a non-negativity constraint on the weights, meaning that the non-zero sparse coefficients should be positive. Is there an easy way to do this? Or another Java package that implements this?


title: Single Multi-label classifier or multiple single-label classifiers?

text: I have a problem to classify my data that can fit into more than one class at same time. Based on an initial study, I came across "Multi-label classifier" that can classify data into more than one class simultaneously. Or as its internal implementation goes, I can have multiple single-label classifiers. I'm fairly new to this domain and have been dealing with single-label classifiers till now. Can someone suggest me which way to opt for based on ease-of-implementation, any reusable existing libraries. Thanks


title: When using the MIXED procedure in SPSS should set my 'time' variable as scale or categorical?

text: I currently have a modelling issue, specifically in determining how to code the 'time' factor in my analysis, when the DV has already been adjusted for age. My DV is a measure of child development, and it has already been adjusted for child age to produce a _z_ -score. Conversion is not linear as the actual assessment varies depending on child age; scores are adjusted based on age of the child on day to the nearest fortnight. Data on this measure were collected at 3 times (1 year, 3 years and 4.5 years of age). This has been coded into SPSS as `Months` and `Visit`. I therefore have 3 measures of 'time' in my data set:   1. `Age` – Child age at assessment (scaled to 2 decimal places),   2. `Months` – General age at which assessment was completed (i.e. 12, 36 and 54 months), and   3. `Visit` – coded as `0`, `1` and `2`. I would like to use the `MIXED` procedure in SPSS to analyze the DV in relation to both 'time' and several other binary variables (e.g. attendance at preschool). I have determined that using the `/REPEATED` command and an AR(1) covariance structure gives the best fit; however, I am not 100% sure how I should code time. Should I use `Age`, `Months` or `Visit`? If I use `Months` or `Visit`, should it be coded as `Scale` (i.e. the `WITH` command) or `Categorical` (`BY` command)? Remember that raw scores have already been adjusted for child age. E.g. 1.               MIXED Z_Score BY Months         /FIXED=Months         /PRINT= SOLUTION TESTCOV         /REPEATED= Months | SUBJECT(UniqueID) COVTYPE(AR1).        E.g. 2.               MIXED Z_Score WITH Months         /FIXED=Months         /PRINT= SOLUTION TESTCOV         /REPEATED= Months | SUBJECT(UniqueID) COVTYPE(AR1).        Thanks for your help.


title: geographic distance and mahalanobis distance

text: I am trying to match individuals based on monthly consumption and geographic consumption in a dense metro area. Essentially I want to create treatment and control pairs with the having both geo location and consumption (total consumption and consumption pattern) as close to the treatment as possible. At first I was just creating a radius of half a mile around each treatment unit and selecting a control that had the closest consumption for winter and summer months—I used a Euclidian distance measure between the treatment and control winter and summer consumption variables. However, this approach doesn’t account for correlation among the variables, so I wanted to measure the mahalanobis distance by first computing the PCA scores for these variables and then finding the Euclidian distance. This got me thinking. Instead of looking at individuals within a particular radius, I could use the standardized longitude and latitudes in the distance calculation. I would still put some kind of constraint on how far away the potential controls can be (I wouldn’t want to compare units form a upper class neighborhood to a ghetto), but the distance metric would include a “closeness” component this time around. Then I thought I could look at each month instead of total winter and summer to better capture total usage as well as consumption patterns. When I did this I found in both the PC1 and PC2 large loadings on the WE- coordinate. Also, I found loadings near zero for winter months and positive for summer months on the first PC and negative loadings on winter months and positive loadings on summer months on the second. I thought this was kind of an interesting pattern. When plotted the monthly consumption of the control and treatment pairs the patterns looked similar, so it seems to have worked. I am just not sure the coordinates have real meaning in this sense. It might be worth estimating the PCA’s on the set of all possible controls instead of in smaller subsets if I am going to used them, but I think some areas are just too different. I am just not sure if there is much sense in including the coordinates. That is, is it incorrect to use longitude and latitude variables in a principal component analysis? I am thinking they should orthogonal variables naturally so the correlation might be misleading. I got similar results using the radius method, but I would like to stick to one method and run with it. Thank you


title: Variance of sample proportion decreases with n but of a count increases with n - why?

text: I've got an intuitive block with this. For a binomial problem, the standard deviation of a count is $\sqrt{np(1-p)}$. Conversely, the standard deviation of the sample proportion decreases with increasing $n$ and is $\sqrt{\frac{p(1-p)}{n}}$. I can do the division by $n$ but I don't have a feel why standard deviations move in opposite directions.


title: understanding of effect of $\alpha$ in Dirichlet distribution

text: When reading the topic modeling tutorial written by Blei, KDD 2011 tutorial I was confused about a set of diagrams which aim to show the effect of $\alpha$ in Dirichlet distribution. For example, for the plot with $\alpha=1$, what am I suppose to discover? What does item mean here? Do those 15 items mean a 15-dimensional probability vector? The $\alpha=1$ is assumed to lead to a uniform distribution. But I feel confused about how to connect this plot with a uniform distribution. I also list other plots with $\alpha=100$ and $\alpha=0.01$. ![enter image description here](http://i.stack.imgur.com/jO3Bo.jpg) ![enter image description here](http://i.stack.imgur.com/QikYy.jpg) ![enter image description here](http://i.stack.imgur.com/BmSoz.jpg)


title: Official name for "symmetric percent difference" function (x-y)/max(x,y)

text: I frequently use this formula to compare two positive numbers $x$ and $y$ to see if they are "more different" than some threshold: $$ x-y \over \max(x,y) $$ It is nice because it is symmetric and bounded to $[-1,1]$ (unlike relative percent difference). I call it a "symmetric percent difference." I see a similar formula on this Wikipedia page, apparently generalized to negative or positive numbers, but it's not named: $$ |x-y| \over \max(|x|,|y|) $$ Does anyone know the official name for this function? Note: Another similar function, bounded to $[0,1]$, is used to calculate sMAPE: $$ |x-y| \over x+y $$


title: In confusion with a Bayesian statistical problem

text: I was learning some probability basics. I am stuck with a problem, that I need your help with in solving. An $e$-fair coin is a coin with probability of head $(\theta)$ in interval $[\frac{1}{2}-e,\frac{1}{2}+e]$. For a Bayesian analysis of the question whether a given coin is $e$-fair, we may assume a uniform prior for $\theta$. Say, we observe $n$ heads out of $n$ tosses of the coin. What is the posterior probability that the coin is $e$-fair? We observe $n$ heads out of $n$ tosses of the coin. What is the minimal number $n$ in order to conclude that the coin is not $0.1$-fair with a probability of at least $99\%$?


title: Differences between intrinsic dimension and correlation dimension

text: I would like to know what is the difference between Intrinsic dimension and correlation dimension. Searching results say the intrinsic dimension = number of variables used to specify the system in euclidean space. So, what is then correlation dimension? I am aware that correlation dimension is a type of fractal dimension. Paper1 talks about intrinsic dimension and Paper2 explains about using maximum likelihood estimation for finding correlation dimension using Takens formula. What is the link between the two dimensions?


title: finding probability density function

text: Suppose $X_1,\ldots,X_n$ are random samples of a $U(0,\theta)$ distribution. If $Y_n$ are the largest order statistics of the sample, and $V=\displaystyle\frac{n\bar{X}}{Y_n}$, how can I calculate the probability density function of $V$?


title: Serial correlation: estimation vs robust SE

text: We have route-level data (that I cannot share) on monthly bus ridership in New York City, creating a panel $N= 185$, $T=36$. We estimate a fixed effects model and random effects model with R's `plm` package. (For this MWE, I use the `Grunfeld` investment data, which illustrates the problems I am seeing fairly well).               library(plm)     data("Grunfeld")     model_1A <- lm(inv ~ value + capital, data= Grunfeld)     fe.plm <- plm(formula(model_1A), model="within", index=c("firm", "year"),                   data=Grunfeld)     re.plm <- plm(formula(model_1A), model="random", index=c("firm", "year"),                   data=Grunfeld)      A Hausman test indicates that the FE model is preferred, because the estimates differ (note that in the MWE, we fail to reject).               phtest(fe.plm, re.plm)          ##  Hausman Test          ##data:  formula(model_1A)     ##chisq = 2.3304, df = 2, p-value = 0.3119     ##alternative hypothesis: one model is inconsistent          pbgtest(fe.plm)          ##  Breusch-Godfrey/Wooldridge test for serial correlation in panel models          ##data:  formula(model_1A)     ##chisq = 65.0632, df = 20, p-value = 1.14e-06     ##alternative hypothesis: serial correlation in idiosyncratic errors          pbgtest(re.plm)          ##  Breusch-Godfrey/Wooldridge test for serial correlation in panel models          ##data:  formula(model_1A)     ##chisq = 69.9495, df = 20, p-value = 1.856e-07     ##alternative hypothesis: serial correlation in idiosyncratic errors      This is where what I think should happen diverges from what many people try to do. My understanding of serial correlation is that it affects the standard errors but not the coefficients. This would suggest to me a serial correlation-robust standard error. For instance (as in this answer),               fe.rse <- sqrt(diag(vcovHC(fe.plm, type="HC1", cluster="group")))     re.rse <- sqrt(diag(vcovHC(re.plm, type="HC1", cluster="group")))      ** Why not just use sc-robust standard errors?** But what many authors do instead is include specific AR(1) or ARMA disturbances because Stata makes this easy. For the FE models, we can use `gls` from the `nlme` package on demeaned data (note: `fe.plm` and `fe.gls` are virtually identical),               # within estimator is demeaned     demean <- numcolwise(function(x) x - mean(x))     Grunfeld.dm <- ddply(Grunfeld, .(firm), demean)     Grunfeld.dm$year <- Grunfeld$year          fe.gls <- gls(update(formula(model_1A), .~.-1), method="ML",  data=Grunfeld.dm)     fear.gls <- update(fe.gls,  correlation = corAR1(form = ~ year | firm))     fearma.gls <- update(fe.gls, correlation = corARMA(form = ~ year | firm,                                                         p=1,q=1))      The RE models can be estimated in with `lme` in the `nlme` package (again, `re.plm` and `re.lme` are identical).               re.lme <- lme(fixed = formula(model_1A), random = ~ 1|firm, data = mta)     rear.lme <- update(re.lme,  correlation = corAR1(form = ~ year | firm))     rearma.lme <- update(re.lme,  correlation = corARMA(form = ~ year | firm,                                                          p=1,q=1))      There are a few things I don't understand about this:   * **Why do the coefficients change when serial correlation doesn't (shouldn't?) affect estimates?**   * **Can we still use the Hausman test to select between FE and RE models with autoregressive errors?**   * **How can we test for residual autocorrelation? And if it exists, wouldn't we _still_ need a robust standard error?**


title: Why maximum likelihood and not expected likelihood?

text: Why is it so common to obtain maximum likelihood estimates of parameters, but you virtually never hear about _expected_ likelihood parameter estimates (i.e., based on the _expected value_ rather than the _mode_ of a likelihood function)? Is this primarily for historical reasons, or for more substantive technical or theoretical reasons? Would there be significant advantages and/or disadvantages to using expected likelihood estimates rather than maximum likelihood estimates? Are there some areas in which expected likelihood estimates _are_ routinely used?


title: Aggregating based on inverse probability weight

text: I have a data set that includes life quality indices for many cities. I want to aggregate the indices to the metropolitan level. Since some cities in each metro area are missing, the data set provides a weighting variable, which was calculated based on the inverse of the probability of a city being in a metropolitan area. For example, in Metro A, 3 cities are in the data set and their indices are 1.2, 0.8, and -0.5. Their weights are 8.2, 3.8, and 6.8 respectively. Here is how I calculated the aggregated index for Metro A: Index for A = 1.2*8.2/(8.2+3.8+6.8) + 0.8*3.8/(8.2+3.8+6.8) + -0.5*6.8/(8.2+3.8+6.8) I wonder if this is the correct way to do it. Any suggestion is appreciated!


title: What is the sample-size range for which the median should be preferred to the mean as a measure of central tendency? And why?

text: I have firing-rate data for neurons ranging from 40-60 in number for various experimental groups. I've seen both mean and median being used in the literature, but without much explanation for the choice. My boss thinks that usually for sample sizes > 4-5, mean is ok to use. But I would like to have a better justification for my choice since the results I get with the two quantities are quite different.


title: Statistical models on Birth interval data

text: Birth interval means time between successive births. I have a birth interval data from DHS survey. If we want to see the effect of familial and socioeconomic co-variates on birth interval. How can we see this? What types of statistical models can be used?


title: dlmForecast function in dlm R package forecasting constant values for seasonal series

text: I have a question regarding the use of the dlm CRAN package for forecasting values of a seasonal time series. I've built a dlm model combining a stochastic local level model with a stochastic trigonometric (Fourier representation) seasonal component of period 96 (measurements every 15 mins with a daily cycle). I used dlmMLE to estimate the parameters for my data and filtered and smoothed the series which all seems to be working fine. However, when I try to use the dlmForecast function to predict out-of-sample observations, the predictions stay constant. The value of all "predictions" are equal to the sum of the filtered level and filtered seasonal components for the final observation in the series. I have used dlmForecast with several other models including a model with a seasonal factor component but never before with a trigonometric seasonal component. I notice in the documentation for dlmForecast it says "Currently, only constant models are allowed" so I wonder if this applies to trigonometric seasonal models.


title: Calculating distances in SAS PROC CLUSTER method=centroid

text: I'm trying to understand the calculations behind the centroid method of cluster analysis in SAS's `PROC CLUSTER`. I am able to calculate the reported distances at each join when I set the options to `method=centroid nonorm`, but I am unable to figure out how the reported distances at each join are calculated when I set the options to `method=centroid nonorm nosquare`. In the documentation, it states that _"Distance data are squared unless you specify the NOSQUARE option."_ Below is some code using a very simple example with 5 points in a plane (from Afifi and Clark, 1996).               data ac;         input x y;         cards;     1 1     1 2     3 6     5 7     8 5     ;          * this gives centroid distances of 1, 2.2361, 4.272, and 6.2472 ;     proc cluster data=ac method=centroid nonorm;         run;          * this gives centroid distances of 1, 2.2361, 3.7933, and 5.0593 ;     proc cluster data=ac method=centroid nonorm nosquare;         run;      Here are my calculations of the distances "by hand" for the first `PROC CLUSTER` statement.               CL4: (1, 1) and (1, 2)         distance = sqrt((1-1)^2 + ((1-2)^2) = 1         centroid = ((1+1)/2, (1+2)/2) = (1, 1.5)          CL3: (3, 6) and (5, 7)         distance = sqrt((3-5)^2 + ((6-7)^2) = 2.236         centroid = ((3+5)/2, (6+7)/2) = (4, 6.5)          CL2: CL3 and (8, 5)         distance = sqrt((4-8)^2 + ((6.5-5)^2) = 4.272         centroid = ((3+5+8)/3, (6+7+5)/3) = (5.3, 6)          CL1: CL4 and CL2         distance = sqrt((1-5.3)^2 + ((1.5-6)^2) = 6.247         centroid = ((1+1+3+5+8)/5, (1+2+6+7+5)/5) = (3.6, 4.2)      Can anyone help me calculate the distances by hand for the second `PROC CLUSTER` statement?


title: Meta-analysis in R with multiple SNPs

text: I'm conducting a meta-analysis on two genome wide association studies (GWAS), each consisting of 150 SNPs, for which I computed summary statistics for association. I'm using the R package meta and the function `metagen`. It seems to me that this function is able to combine studies but on single measurement effect (I have 150), one for each study. For example, I may have 3 SNPs on each study               STUDY    SNP    OR    SE     GWAS_A  rs694739    0.6691  0.07588     GWAS_A  rs9858968   0.1091  0.01588     GWAS_A  rs1529267   0.9291  0.02588     GWAS_B  rs694739    0.6128  0.37344     GWAS_B  rs9858968   0.0332  0.27344     GWAS_B  rs1529267   0.3481  0.81284      so I would like an output in a table format (e.g., as in PLINK), with 3 rows and the corresponding pooling statistics. I know I could manage it with `for` loops or by grouping with the `by` function on the summary output of `metagen`, but I wonder whether there is another R package/function to easily accomplish this.


title: Attainable correlations for exponential random variables

text: What is the range of attainable correlations for the pair of exponentially distributed random variables $X_1 \sim {\rm Exp}(\lambda_1)$ and $X_2 \sim {\rm Exp}(\lambda_2)$, where $\lambda_1, \lambda_2 > 0$ are the rate parameters?


title: Pregibon test for linearity vs. Ramsey's RESET test

text: Does every Ordinary Least Squares (OLS) regression model have to pass both the Pregibon Test for Linearity (sometimes called link-test) and the Ramsey RESET tests? I am working on an OLS model and it passes the Pregibon test when I include all the variables I wanted without squared terms. However, it fails to pass the RESET test. Do I have to include quadratic terms or is the Pregibon test "enough" to claim that the model is well specified? Ideally would prefer to keep the model simple, so to not include quadratic terms. thanks!


title: Maximum value of Chi squared while keeping marginal frequencies

text: Let's say we have a RxC contingency table, with R=C. How could one best repopulate the table such that:   1. all counts reside in the off-diagonal cells (zero's on the diagonal);   2. marginal frequencies remain the same;   3. the value of Chi square is maximal. I can imagine a very 'exhaustive' approach by calculating the Chi square for each possible repopulation and taking the table that produces the highest Chi square value, but there has to be an elegant much better way to solve this! A formal proof may be difficult, but a smart algorithm should exist.


title: Using simpler models in place of more generalized and complex models

text: I was reading about BATS (Box-Cox transformation, ARMA errors, Trend and Seasonality) and TBATS (Trigonometric, Box-Cox transformation, ARMA errors, Trend and Seasonality) models. I was wondering whether it's possible to get better forecasts with simpler models like standard Holt Winters than using BATS or TBATS?


title: Interpreting negative beta coefficients

text: I'm running a negative binomial regression to test the intensity of violent events in a given location. One of my X variables has a negative beta coefficient (BC), however the hypothesis behind that variable is a 'negative decline in X contributes to a positive chabe in voilence'. My X variable is vegetation cover - so, if there is a decline in vegetation, the hypothesis would posit an increase in violence. I was wondering, because the hypothesis behind that variable is for a negative, rather than positive change in X effecting a postive change in Y, would I interpret a negative BC as supporting the hypotheis or not? Am I right in thinking the regression would quantify the BC through looking at a postive change from the mean? I don't know how to interpret this, unless there is a way for me to change my varibale/ if I need too?


title: Best Practices for Creating 'Tidy Data'

text: Hadley Wickham wrote a stellar article called "Tidy Data" (link) in JSS last year about data manipulation and getting the data into an "optimal" condition in order to perform analysis. However, I was wondering what were best practices in terms of presenting tabular data in a work setting? Let's say your coworker asks you to provide him with some data. What are some general rules you use when structuring that data? Are the guidelines in "Tidy Data" just as applicable in instances where you're sharing data with non-data professionals? Obviously, this is very context-specific but I'm asking about high level 'best practices'.


title: How to replace only SOME missing data

text: I have two types of missing data (non-responses and "Don't know or does not apply") and I want to replace only the non-responses (with mean or by using EM or by using Multiple Imputation) and to see if this is a good solution. How can this be done in SPSS? Even in the simplest case, when replacing missing values with mean for example, I can't find a way to impute only some of the missing data, in this case the non-responses?


title: How to evaluate a statistical test for fairness of a coin?

text: I suspect that a given coin is not fair, i.e that the probability of heads is not $\frac{1}{2}$. My question is not how we can prove the unfairness, but if we derive a test to check the fairness, How we can evaluate the test? Say, for $n$ number of tosses, the coin always turned up heads. What is the minimal number $n$ in order to reject the fairness hypothesis at significance level of 1%? Thank you in advance for your attention.


title: Bayes specify informative and uninformative prior

text: I would be really thankful for any hint Let $x_1,\dots x_n$ be a sample from a geometric distribution with probability parameter $\theta$ $f(x_i\mid \theta)=\theta(1-\theta)^{x_i}$ The likelihood is given by $L(\theta\mid x_i, \dots, x_n)=\prod_i f(x_i\mid \theta)$ Assume we have observed a sample of size n=5 with values $x_i=5,6,9,5,4$. Specify a) an informative prior, (b) and uninformative prior, and determine the posterior distribution for each setting. Make a plot of prior, likelihood and posterior distribution.


title: Equations of significance probabilities

text: Consider a population of independent light bulbs with an exponential lifetime distribution with mean $\mu$. It is claimed that their expected lifetime is 1000 hours. A definition of a **100(1−)%** confidence interval obtained from an observation t0 is the set of all $\mu_0$ which are not rejected in a test of a null hypothesis $\mu$= $\mu_0$ against an alternative hypothesis $\mu$$\neq$$\mu_0$. * * * One particular light bulb fails after 622 hours. Solve the equations of the two significance probabilities Pr(T ≥ 622 |$\mu_0$) = 0.05 (for a test of $\mu$= $\mu_0$ versus $\mu$> $\mu_0$) and Pr(T ≤ 622 | 0) = 0.05 (for a test of $\mu$= $\mu_0$ versus $\mu$< $\mu_0$) for $\mu$. Determine the range of values of such that both of the probabilities Pr(T ≥ 622 | $\mu$) and Pr(T ≤ 622 |$\mu$) are at least 0.05. (This range gives an equi-tailed 90% confidence interval for $\mu$.) * * * What inference may be drawn from this interval about the claim that the expected lifetime is 1000 hours?


title: Why is Bayesian error concave down as a function of prior P(w1)?

text: Given a two categorization problem, why Bayesian error is concave down as a function of prior $P(w_1)$? $w_1$, $w_2$ are two underlying concepts with prior $P(w_1)$ and $P(w_2)$. We observe the evidence $x$, which is related to $w_1$ and $w_2$ through conditional probability $p(x|w_1)$ and $p(x|w_2)$. Based on $x$, we will do action $a_1$ or $a_2$, and the cost according to each action is $c(a_i|w_j)$, i.e. $c(a_1|w_1)=0$ (guess correctly) and $c(a_2|w_1)=1$ (guess wrong). And finally Bayes error is the expectation of $\text{cost}(\text{action}|x)$.


title: Ordered combinations

text: I have a combinatorics question. Say you have two sequences: $$X_{1},X_{2},X_{3},\ldots,X_{N_{1}}$$ and $$Y_{1},Y_{2},Y_{3},\ldots,Y_{N_{2}}$$ How can i pair up elements from the $X$'s with elements of the $Y$'s, such that if I pair up $Y_{1}$ with $X_{3}$, then I can only pair up $Y_{2}$ with $X_{n}$ such that $n$ is greater than the index of the $X$ paired up with the previous $Y$, in this case $n>3$.


title: Discretization of skewed data (time durations)

text: I have data that describes the duration of how long a person views a webpage. This is quite varied and in the context wherein I gathered the data, it was very skewed. People mostly spent short amounts of time in a webpage but sometimes spent a significant amount of time viewing it. I want to discretize the durations (in seconds) into short, medium and long but I don't know how I should do this if the data is skewed. Initially, I just used tertiles but it seemed kind of off. Tertiles assume equal membership but I'm not sure if this is right because of the skew. Any ideas on a better way of discretizing the values?


title: Calculating card deck shuffling odds

text: Start with a standard deck of 52 cards in some known order; each card can be tagged with its ordinal position in the deck. Then apply a random shuffling algorithm (i.e. Knuth's shuffle). After the shuffle, it is possible for one or more cards to end up in the same ordinal position where they were before the shuffle -- i.e. the shuffle did not move them. How can I calculate the probability of this occurring based on the number of unmoved cards? In other words, calculate _P(n)_ as the probability that, after the shuffle, at least _n_ cards have not moved.


title: Why is k called representer of evaluation in the definition of kernel functions

text: Why is $k$ called _representor of evaluation_? From the book "Learning with kernels" by Bernhard Schölkopf we have the following lines (page 33): > $\langle k(.,x),f\rangle = f(x)$, in particular $\langle k(.,x), > k(.,x')\rangle = k(x,x')$ According to the book this interesting property of $\phi$ follows from definition. How? I am unable to understand this and this is crucial for understanding the concept of _reproducing kernel Hilbert spaces_. Any help appreciated. Luckily, the part of the book which needs to be referred _i.e._ section 2.2.2 (starts on page 32) is a part of the preview in Google books. Also note that this section is independent of other sections.


title: Why do we divide?

text: This may seem like a rather simple concept, but I'm trying to get an intuitive understanding of why we divide. Common formulas in statistics that use division are the calculation of the t statistic and variance. A less well known calculation in evolutionary biology that uses division is strength of female preference for males: $$\frac{\text{Detectability} \times \text{Honesty} × \text{Benefit from male}}{\text{Costs of exerting choice}}$$ So could someone provide a general, high-level explanation of why we use division?


title: There are two students born exactly in one day of year. How to prove?

text: I remember once in our Statistics class, our professor mentioned an interesting problem, that is: > In our class, consisting of 30 students, I bet that there are two students > born on an exact day of year. It was quite clear to all of us that he will lose the bet, but then he asked all the students to write out their birth date on the board and surprisingly three of us were born in the same day of year! It was very interesting to us. He wrote the proof to his statement, showing that If there are 40 students in one classroom, the probability that at least two students are born in the same day of year is over 90 percent (as I recall). A few days ago, I told my friends about this, they all laughed at me saying it is impossible. I tried to give them the proof but I failed. I've forgotten everything! I guess it is a famous problem in statistics, can anyone provide the proof or any link to any docs providing the proof to this statement?


title: Bonferonni Correction interaction

text: Suppose we have two variables $x_1$ and $x_2$ and an interaction term $x_1 \cdot x_2$. Suppose we set the family-wide error rate to $\alpha = 0.05$. For the Bonferonni correction, would be look at $\alpha/2$ or $\alpha/3$?


title: Hessian matrix and initial guess in logistic regression

text: The log-likelihood function for logistic function is $$l(\theta) = \sum_{i=1}^m(y^{(i)}\log h(x^{(i)}) + (1-y^{(i)})\log(1 - h(x^{(i)})))$$, where $$h(x^{(i)}) = \frac{1}{1 + e^{-\theta^Tx^{(i)}}}\,.$$ In order to obtain maximum likelihood estimation, I implemented fitting the logistic regression model using **Newton's method**. I encountered 2 problems:   1. I try to fit the model to my data, but during the iterations, a _singular_ Hessian matrix is encountered, what do I do with this kind of problem?   2. With different initial guess $\theta$, will the model converge to different results?


title: General-to-specific or ordered response probit?

text: I have a database with the Chilean equivalent SAT for 2010. It shows the highest score achieved by one student that belongs to one school. The scores goes from 100 to 820, and it has been published as a 1-to-100 school ranking. Database also contains several qualitative information of schools (managerial issues, building facilities, staff experience, etc.) If I want to obtain the main determinants inherent of one school that allows students to be at top ranking, do I use general-to-specific estimation, an ordered response probit, or other?


title: Expressing ridge solutions as sum of vectors

text: I'm having trouble following the derivation of the ridge solutions as a sum of vectors in formula 3.47 of The Elements of Statistical Learning on pg. 66. It shows the following: > > $\mathbf{X}\hat{\beta}^{ridge} = \mathbf{X}(\mathbf{X}^T\mathbf{X} + > \lambda\mathbf{I})^{-1}\mathbf{X}^T\mathbf{y}$ >> >> $\mathbf{X}\hat{\beta}^{ridge} = \mathbf{U}\mathbf{D}(\mathbf{D}^2 + \lambda\mathbf{I})^{-1}\mathbf{D}\mathbf{U}^T\mathbf{y}$ >> >> $\mathbf{X}\hat{\beta}^{ridge} = \sum{\mathbf{u_j}\frac{d_j^2}{d_j^2 + \lambda}\mathbf{u_j}^T\mathbf{y}}$ >> >> where $\mathbf{X} = \mathbf{UDV}^T$ I'm having trouble seeing how to get to line 2 after simplifying line 1 with $\mathbf{X}\hat{\beta}^{ridge} = \mathbf{UDV}^T(\mathbf{VDU}^T\mathbf{UDV} + \lambda\mathbf{I})^{-1}\mathbf{VDU}^T\mathbf{y}$ and then: $\mathbf{X}\hat{\beta}^{ridge} = \mathbf{UDV}^T(\mathbf{VD ^2V} + \lambda\mathbf{I})^{-1}\mathbf{VDU}^T\mathbf{y}$ because $\mathbf{U}^T \mathbf{U} = \mathbf{I}$ Also, is there a general theorem that helps yield line 3? I can't wrap my head around going from matrix products to vector sums. It's been a decade since I've taken linear algebra, so if you have any good resources for dealing with this type of math, please let me know.


title: How can I explain these linear regression charts well on my scientific poster?

text: I know that these laboratory analysis reports each have a linear regression relationship and two have a positive slope and one has a negative slope. I am taking my first statistics class and want to be able to explain these very well in a research poster I am doing for work. Am I missing any important information? The negative slope indicates that the values are decreasing together negatively. The positive means they are increasing together. The ADL concentration chart does not have a very good relation but there is still a relation. ![ADL linear regression](http://i.stack.imgur.com/xxZYO.png) ![ADF linear regression](http://i.stack.imgur.com/F1mSM.png) ![r^2, p-value](http://i.stack.imgur.com/vWwFv.png) ![NDF linear regression](http://i.stack.imgur.com/qxpuk.png)


title: Should EFA be repeated after deletion of items with low factor loadings?

text: After running Exploratory Factor Analysis in spss v.17 I deleted 4 measurement items which yielded low factor loadings. Now I want to proceed to Confirmatory Factor Analysis to calculate AVE and CR. My Question is should I repeat the EFA to obtain the new factor loadings after deletion of those 4 variables(items) and conduct the CFA based on them or i should do it with the original factor loadings? Thnx 4 your help.


title: Intuition on learning rate or step-size for perceptron algorithm

text: Recall the perceptron algorithm: _cycle through all points until convergence_ $\text{if }\, y^{(t)} \neq \theta^{T}x^{(t)} + \theta_0\,\{\\ \quad \theta^{(k+1)} = \theta^{k} + y^{(t)}x^{(t)}\\ \}$ I was studying a modification to to the update rule such that the new update rule is: $\theta^{(k+1)} = \theta^{k} + \eta_k y^{(t)}x^{(t)}\\$ where: $\eta_k = \frac{ Loss (y^{(k)} \theta^{(k)} \cdot x^{(k)} ) }{\left \| x^{(k)} \right \|^2}$ and the loss function was the hinge loss. i.e: $Loss(y^{(k)} \theta^{(k)} \cdot x^{(k)}) = max\{0, 1-y^{(k)} \theta^{(k)} \cdot x^{(k)}\}$ I was trying to understand the new weight $\eta_k$ and understand why it was the way it was. I think intuitively I can see why the hinge loss is being used because the least confident we are about our prediction, the higher the loss value it will give it (since we are more concerned of correcting that specific example I guess...), however, I was not sure what the denominator was doing. It seems to me it's some kind of an attempt to normalize the step-size or the weight, but was unsure how to interpret it. However, I was not 100% sure why the numerator was the way it was and any additional insight on either/both would be appreciated! Thanks!


title: Is it OK to fit random effects before fixed effects?

text: I'm using lmer to analyse my data, building nested models and using anova() to compare them against each other in an incremental way. Now, I know enough to only test a single term at a time (i.e. only one term changes between model 1 and model 2, so that when I compare them with anova I know what question I am answering), my question rather relates to the order in which you should test your terms of interest. **Is there a preference or a rule that states whether you should fit your fixed effect terms before your random effect terms? If so, what is the reason?** Thanks all very much


title: Negative Binomial "Process"

text: I wish to model the number of bugs caused by software development. This is intuitively sort of a Poisson process, however it is overdispersed. One thing we can do in this case is to use a negative binomial distribution (because negative binomial approaches Poisson as `r` gets larger, or because we might think the parameter $\lambda$ of the poisson is itself gamma-distributed.) I'm not sure how to do this though. For example, we have that $$\lim_{r\to\infty}NB\left(r,\frac{\lambda}{\lambda+r}\right)=\text{Poisson}(\lambda)$$ Given that a poisson process of duration $t$ can be modeled as $\text{Poisson}(\lambda t)$ I guess we could look at $NB\left(r,\frac{\lambda t}{\lambda t+r}\right)$ - is that correct? Given that I know $t$, it seems like I should be setting $r=t$. At a more technical level, `glm.nb` from the `MASS` package seems to fit $r$ not the dispersion parameter and I don't see an obvious parameter to change this. Any insight at the theoretical or technical level would be appreciated.


title: Is it possible to visualize my data set in Excel?

text: My data set consists of columns containing universities in the Boston area, and rows containing zip codes. For each university and zip code, there are two data points: number of students in campus housing, and number in non-campus housing. So that's two independent variables `(x, y)`, and for each such pair, there are two dependent variables `z1` and `z2`. I'm trying to visualize this in Excel, but not finding a chart type that could handle it. My best idea would be a sort of bubble chart, with x and y on the axes, and the z1 and z2 represented as bubbles, but I can't make this work. Can anyone suggest a way to do this, or another visualization that would work?


title: Chi-squared test when two vectors have different lengths

text: I want to compare two groups (vectors) having different lengths to see if their distributions differ considerably. However, when I use `chisq.test`, I get an error due to inequalities of lengths:               > chisq.test(factor(Group1), factor(Group2))      results in:               Error in chisq.test(factor(Group1[, col]), factor(Group2[, col])) :      'x' and 'y' must have the same length      Could someone help me with this? Thanks.


title: Interpreting standardized betas with log-transformed dependent variables

text: I've run a multiple linear regression using standardized Betas and a log- transformed dependent variable; the latter transformation is primarily for better approximating linear regression assumptions but it also produces a nifty interpretation (described below). Standardized Betas provide the interpretation: a 1-SD increase in the predictor is associated with a Beta X SD change in the outcome. Log-transformation of the outcome variable allows for the interpretation: a 1 unit increase in the predictor is associated with a Beta times 100% change in the outcome. Are the Betas from my log-transformed and beta-standardized model interpreted as: a 1-SD increase in the predictor is associated with a Beta times 100% increase in the SD of the outcome? Or is there a mathematical nuance that changes this all-together?


title: Measuring correlation when multiple data points in a sample are provided by the same subject

text: I want to measure the correlations between ratings on different items of a scale. In my sample, there are multiple data points from the same subjects, for each subject to the same extent. I am able to trace the ratings back to the subject who made them. What is the appropriate correlation measure accounting for these dependencies? I might need a non-parametric one. Here's a toy example: Subjects have to rate the acidity (item 1) and consistency (item 2) of fruit. Thus, they are shown different stimuli (apple, oranges, pears, ...) and rate for each fruit once on each item. Now I want to examine whether the two items acidity and consistency correlate somehow. Nevertheless, If I want to embrace all stimuli in my observation, there will be multiple ratings from the same subject. Nevertheless, there is only one rating per item, stimulus and subject.


title: Non-parametric statistical test with one measurement variable and two nominal variables

text: I have the following data: feces density of a carnivore species subdivided according to sex. I am looking for appropriate test suggestions. I have calculated the density according to main habitat types, topography, and road types. To determine if results are statistically significant, I would believe the Friedman`s test or Wilcoxon signed-rank test would apply in this case since there is one measurement variable (feces density) and two nominal variables (e.g., habitat type and sex; road type and sex; topography and sex) and both are non-parametric statistical tests. Any help would be greatly appreciated. Thank you in advance!


title: Use hierarchical clustering in R to cluster items into fixed size clusters

text: I am trying to use R to do Kmeans clustering and as most people I ran into the challenge of determining when to finish. I have 10,000 items and potentially 10 times of that down the road. My goal is to create a series of clusters with minimal size (e.g. 50 items per cluster) OR reasonably similar items. In other words, I don't want any of my output clusters to be too small (even if the items are quite different from each other), but I also don't mind if the clusters are too big as long as the items are similar enough. I imagine I can use some kind of divisive hierarchical approach. I can start by building a small number of clusters and examine each cluster to determine if it needs to be split into more clusters. I can keep doing this till all clusters meet my stopping criteria. I wonder if anyone knows good information on how other people do this?


title: The junction tree theorem

text: In [1] (Page 31,equation 2.12 ) it is claiming that in a graph which is processed by the junction tree algorithm, the joint distribution of the variables could be found by $$ p(x_1, ..., x_m) = \frac{ \prod_{C \in \mathcal{C}} \mu_{C}(x_C) }{ \prod_{S \in \mathcal{S}} [\mu_{S}(x_S)]^{d(S)-1} } $$ In which $d(S)$ denote the number of maximal cliques to which it is adjacent. The problem is that, I don't really see where this distribution is coming from. Any ideas? [1] http://www.eecs.berkeley.edu/~wainwrig/Papers/WaiJor08_FTML.pdf


title: Finding the right model

text: I have a data frame with several predictors, lets call them `pred1` through `pred3`, and a result column. Now I need to specify the right model. I could randomly try:               svm.model <- svm(result ~ pred1+pred2+pred3,       data = train)     # or     svm.model <- svm(result ~ pred1*pred2+pred3^2,     data = train)     # or     svm.model <- svm(result ~ log2(pred1)+pred2*pred3, data = train)     # etc. etc.      But there must be a better automatic approach in R to model selection?!


title: What is the loss function of hard margin SVM?

text: People says soft margin SVM use hinge loss function: $\max(0,1-y_i(w^\intercal x_i+b))$. However, the actual objective function that soft margin SVM tries to minimize is $$ \frac{1}{2}\|w\|^2+C\sum_i\max(0,1-y_i(w^\intercal x_i+b)) $$ Some authors call the $\|w\|^2$ term regularizer and the $\max(0,1-y_i(w^\intercal x_i+b))$ term loss function. However, for hard margin SVM, the whole objective function is just $$ \frac{1}{2}\|w\|^2 $$ Does that mean hard margin SVM only minimize a regularizer without any loss function? That sounds very strange. Well, if $\frac{1}{2}\|w\|^2$ is the loss function in this case, can we call it quadratic loss function? If so, why the loss function of hard margin SVM becomes regularizer in soft margin SVM and make a change from quadratic loss to hinge loss?


title: Comparison of two curves

text: I have two measurements of something. You can think of them both as their own curve with known error on each point. They look quite similar, and if I calculate the R2 value it comes out to >0.9. But what I want to be able to calculate is a P value comparing these two curves (i.e., what is the probability that the difference i'm looking at is just due to noise?). Now I could easily do a student t-test at each point and come up with a P-value at each point, sure. But is there some way to come up with an over-all P value that uses all the points and not just one? Thanks very much for any help.


title: Estimating workers within Industry Sector using commuter data (combining probabilities?)

text: I have some detailed data about the commuting patterns of workers. If I know the following 3 things...   1. The total number of workers commuting from one block to another.   2. The total # and % of manufacturing-sector workers in the origin block.   3. The total # and % of manufacturing-sector workers in the destination block. ... Is it possible to estimate the number of _manufacturing-sector workers_ who commute from the _origin block_ to the _destination block_? Data structure is below. I think it's a probability combination, but I'm too rusty to come up with an elegant solution. Thanks a bunch for taking a look. If successful, this analysis will be used to characterize the regional commuting patterns of specific industry clusters. * * * ![Example structure](http://i.stack.imgur.com/4rAXs.png)


title: Should I control for random effects of participant in an individual differences design?

text: I'm trying to analyse a survey study in which I'm interested in the way that individual differences between my participants influence how they respond to my stimuli. The stimuli are pieces of writing randomly selected from a large pool of such, and I'm not currently interested in any fixed effects of these stimuli. My participants each view 5 pieces of writing and make a binomial choice about each on a single-item DV. I have personality variables from each of the participants that I'm interested in as a predictor of their response. When setting up a mixed-effects model in lmer (with a long dataset) I know that I should include a term for the random effects of the stimuli. My question is whether I should also include a term for the random effects of participant, along with the fixed effects of the personality variable. Because differences between participants on my personality measure are my variable of interest, does it make sense to control for random participant effects? Also, if I do control for the random effects of participant, what am I conceptually doing here? Partitioning out all the participant-level variance that can't be explained by my personality variable? Controlling for the random effects of participant seems to consistently reduce the predictive power of any personality variable I include in the model, which is why I ask. Some illustrative lmer syntax for the two options: Without participant random effects: `glmer(Choice ~ Personality + (1| StimID), data=choice.long, family='binomial')` With participant random effects: `glmer(Choice ~ Personality + (1| StimID) + (1| PartID), data=choice.long, family='binomial')` Am I specifying the latter model right?


title: Is there a relationship between the median of a function of random variables and the function of the median of random variables?

text: **Background** notation: RV= random variable, $\mu=$ mean $m=$ median Jensen's Inequality considers the relationship between the mean of a function of an RV and the function of the mean of an RV. If $f(x)$ strictly convex: $$\mu (f(x)) > f(\mu (x))\mathrm{\hspace{20mm}(1)}$$ Conversely if -f(x) is strictly convex: $$\mu (f(x)) < f(\mu (x))$$ An analogous property of the median has been presented (Merkle et al 2005, pdf). **motivation** I have a nonlinear function of positive random variables. In practice, I find that the function of the medians provides a much better estimate of the median of the function than does the estimate of the mean of the function from the function of the means. I am interested in learning the conditions for which this is true. **question** Under what conditions will the function of a median be closer to the median of a function than the mean of a function is to a function of the mean? Specifically for what types of $f(x)$ and $x$ is $$\mu (f(x)) - f(\mu (x)) > m (f(x)) - f(m (x))$$ **simulation results** I used an empirical approach (the one I know) to investigate this question for a function of a single variable: Interestingly, for $x>0$, $$m(x^2)\simeq m(x)^2$$               set.seed(1)     x<-cbind(rlnorm(100, 1), rbeta(100, 1, 5), rgamma(100,0.5,0.5))     quad <- function(x)x^2          median.x <- apply(x,2,quantile,0.5)     mean.x <- apply(x,2,mean)          colMeans(quad(x))     quad(mean.x)          apply(quad(x), 2, quantile, 0.5)     quad(median.x)      For a slightly more complicated function, my proposal (equation 1) is true               miscfn <- function(x) 1 + x + x^log(x^2) - exp(-2(x)*5^x           colMeans(miscfn(x))     miscfn(mean.x)          apply(miscfn(x), 2, quantile, 0.5)     miscfn(median.x)          abs(apply(miscfn(x),2,mean)-miscfn(mean.x)) > abs(apply(miscfn(x), 2, quantile, 0.5) - miscfn(median.x))      However, before I begin to use this observation in my work, I would like to know more about its conditions. **References** Merkle et al 2005 Jensen's inequality for medians. Statistics & Probability Letters, Volume 71, Issue 3, 1 March 2005, Pages 277-281


title: MCMC Standard error

text: I need to make a simulation to calculate the integral of a pdf $f(x|D)$ over a region $T$.That is, I need to evaluate $$\int_T f(x|D)\; dx$$ to do this, I just defined the function $g(x)=1$ if $x \in T$ and $g(x)=0$ if $x$ is not in $T$. So the integral that I want to evaluate is $$\int_{D(f)} g(x)f(x|D)\;dx$$ where $D(f)$ is the domain of $f$, right? The Monte Carlo estimative of my integral is given by $$\bar{I}=\frac{1}{N} \sum_{i=1}^N g(x_i), x_i \sim f$$ right? And the Monte Carlo Standard Error is $$SE= \sqrt{\frac{1}{N}\sum_{i=1}^N(g(x_i)-\bar{I})^2}$$ I made a Python function to calculate this:               def MCSE(g, estimativa):          #Size of the vector g (which is composed of 0 and 1)     N = int(len(g))     #Auxiliar vector vet     vet = []     for i in range(N):         vet.append(float((g[i]-estimativa)**2))     #Calculate the standard error and return     SE = sqrt(np.mean(vet))     return SE      Is my function correct? Because when I estimate the integral and it gives a value close to the real value of the integral (which I know from other ways), it is giving values close to 0.5, which doesn't seems correct. Thanks in advance!


title: comparing ARMA coefficients

text: I am comparing two sections of a single time series on the date that a policy change has occurred. How do I compare two ARMA coefficients? Do I use t test? Thanks a lot.


title: How to calculate variance - covariance matrix of a matrix?

text: For example, we have an NxP matrix with N rows and P variables. Then we need to calculate a PxP sample covariance matrix. How do I do that?


title: Gradient descent oscillating a lot. Have I chosen my step direction incorrectly?

text: I'm trying to run a basic gradient descent algorithm with a absolute loss function. I can get it to converge to a good solution by it requires a much lower step size and more iterations than had I used square loss. Is this normal? Should I expect absolute loss to take a longer time to come to a good solution or potentially oscillate around a solution more than say squared loss?


title: How to optimize neural words embedding (Colbert & Weston, 2008) vectors?

text: AFAIK, the model works in a following way: We have $S$ - all the sentences in our unlabeled copra.   1. choose $N$ (~100000) most popular words from our unlabeled data = $PW$   2. chose the hyperparameters $d_{dim} \in [50...200]$ - the number of dimensions in our embedding and $w_{size}$ - window size   3. initialize the $LT \in M^{d_{dim} \times N}$ (lookup table) matrix with the random numbers - it associates every word with it's $d_{dim}$-dimensional embedding.   4. define the score function (NN in fact) $score(sentence) = U^T f(W \times concat([LT(w_i) | w_i \in sentence]) + b) $, where $w_i \in sentence$ - every single word in the sentence, $concat$ - is the function that concatenates all vectors given, $f()$ - some nonlinear function (e.g. tanh).   5. optimize function : $ L = \sum \limits_{sentence \in S} \sum \limits_{s \in wnds(w_{size}, sentence)} \sum \limits_{word \in PW} max(0, 1-score(s) + score(s^{word}) $ where $wnds(w, s)$ - all windows of size w in s and $s^{word}$ - corrupted sentence - with central word replaced with $word$. So we maximize the score of "good" (real) sentence pieces over the "bad" (corrupted) ones. The way of optimizing $W, b, U$ - is trivial, just like any other network: backpropagation to whatever - it is written in every paper. But the question is: in fact we must optimize the $LW$ to get embedding (the LW itself), isn't it? So what should we do to find it? Optimize all those $N$ embedding vectors (${d_{dim} \times N}$ numbers) together with $W, b$ and $U$ ?


title: existence of fractional moment of cauchy distribution

text: For the standard Cauchy distribution, do fractional moments exist? If $Y \sim C(1,0)$, is it possible to evaluate $E(Y^{1/3})$?


title: Relationship between Binomial and Beta distributions

text: I'm more of a programmer than a statistician, so I hope this question isn't too naive. It happens in sampling program executions at random times. If I take N=10 random-time samples of the program's state, I could see function Foo being executed on, for example, I=3 of those samples. I'm interested in what that tells me about the actual fraction of time F that Foo is in execution. I understand that I is binomially distributed with mean F*N. I also know that, given I and N, F follows a beta distribution. In fact I've verified by program the relationship between those two distributions, which is               cdfBeta(I, N-I+1, F) + cdfBinomial(N, F, I-1) = 1      The problem is I don't have an intuitive feel for the relationship. I can't "picture" why it works. EDIT: All the answers were challenging, especially @whuber's, which I still need to grok, but bringing in order statistics was very helpful. Nevertheless I've realized I should have asked a more basic question: Given I and N, what is the distribution for F? Everyone has pointed out that it's Beta, which I knew. I finally figured out from Wikipedia (Conjugate prior) that it appears to be `Beta(I+1, N-I+1)`. After exploring it with a program, it appears to be the right answer. So, I would like to know if I'm wrong. And, I'm still confused about the relationship between the two cdfs shown above, why they sum to 1, and if they even have anything to do with what I really wanted to know.


title: Comparison of the tails of two sample distributions

text: I have two set of data that are roughly centered around zero but I suspect that they have different tails. I know a few tests to compare the distribution to a normal distribution, but I would like to compare directly the two distributions. Is there a **simple test to compare the fatness of tail of 2 distributions**? Thanks   fRed


title: Why does maximum likelihood estimation not work in estimating signal in deterministic chaotic noise

text: I have few conceptual questions related to application of chaos in communications. In few application such as radar Chaotic signal reconstruction with application to noise radar system, cryptography, the noise considered is chaotic in nature. Why does MLE not work?


title: How to find Coefficient of Variation with E(X), E(X^2) and sample size?

text: Given E(X), E(X^2) and sample size. How do I find Coefficient of Variation? Variance = E(X^2) - E(X)^2 sd = square root of variance Coefficient of Variation = sd / E(X) however, with the values that were given, I can't get the answer I want. I wonder if the above formula is correct.


title: An approach to visualise/investigate multi-dimensional data

text: I'm attempting to analyse a large chunk of empiric measurements subject to two parametrised transformations. In essence, the functions take 3 'count' parameters - and returns a sequence of floats in each configuration. I'm expecting (hoping) to see some interesting patterns emerge when appropriate parameters are selected. I anticipate that the patterns might be relative between sequences returned for each function - and/or relate to patterns of some kind in the parameters. In case it's relevant, the 3 'count' parameters roughly correspond to:   * 'window size' on the underlying data over which summary statistics are calculated   * A number of consecutive windows used to compute a single summary statistic (i.e. the trade-off between greater spatial or greater temporal accuracy)   * A 'minimum age' - an offset into history of the underlying data. The summary statistics are non-trivial but are independently sensitive to all three parameters. I'm interested in visualisation techniques - suited to ad-hoc enquiry - that will help me experiment with this multi-dimensional data. **EDIT** (In response to Peter Ellis): Your 'pun' is, essentially, my problem. I've not stored the data - instead, I can calculate it 'on the fly' from (otherwise opaque) bulk empirical data. If I were to store the results of my calculation in matrix form, it would be a 4D matrix... The first three dimensions are all temporal - as identified above - then the index for the resulting sequence is also temporal - and the values are a dimensionless problem-specific scalar metric... represented using a floating point number. One might imagine the data to visualise as associating a distribution (function) with each 'voxel' in a cuboid. I visualise these sequences (representing the distribution associated with each voxel) in isolation, as regularly sampled continuous functions. One way I've imagined the data is with the sequences as line-graphs, and three 'tweaking knobs' I can use to 'tune in' a result... another is as any one of thee surface plots with two 'tweaking knobs'. A more elaborate visualisation animates a range of values for one of the 'tweaking knobs' - leaving me one parameter to set manually. The structure of the data does not suggest which parameters are best in which context - and the answer to this may depend upon both the empiric data (yet to be collected in full) and the scale of the parameters. A significant complication is that I don't - as yet - know how sensitive the results of my calculations will be to the parameters above... One objective is to find values for the first two parameters to minimise the impact of the third. Beyond that, and analogies between my parameters and time, I'm afraid the data is pretty abstract... there's no real-world object from which to draw inspiration.


title: Solution for long-term electricity prices forecasting

text: I'm trying to adopt a solution for long-term electricity annual prices forecasting (depending on past electricity prices, past oil prices, past consumption data, etc.) I'm considering some solutions:   * off-the-shelf specific software: Alyuda, Aleasoft, ...   * already-built modules of known software (Excel, MatLab, R): neuroXL for Excel, modules for Simulink, ...   * code from scratch (Python, R) by using known models and Machine Learning: scikit-learn, Weka, regression, ... I would like to find more solutions to test and, if possible, some experiences using them. Thank you very much in advance.