title: CI for conditional probabilities

text: I am searching for some information about confidence intervals for conditional probabilities I thought about treating them as normal probabilities with the formula $$P(Y=i|x=t) \pm z \left(\frac{ \mathrm{Var}(P(Y=i))}{n} \right)^{1/2}$$ with $n$ the sample size of observations $t$. But is this valid?


title: Machine learning on big data: capability of generalization

text: I have being working applying different ML Algorithms oriented to Big Data and I have some open questions that I find interesting to think about. One of the first lectures about statistics begins with the sentence "Since the whole population data is not available, we must take representative samples so that our conclusions on this sample could be extended to the whole population." And Statistics is all about this. However, Big Data allows not only access to the whole population but also the possibility of processing that information (something that was assumed false in statistics since the beginning). So here are my questions:   * If we can process all the data (the whole data), what's the point of doing nowadays cross-validation, sampling and the like?   * If the answer to the latter question is 'overfitting', what if my algorithm is 'overfitted' with millions and millions and millions of data? In my opinion, if my algorithm is able to generalize that much I don't care if it is overfitted.   * Do you think there are limitations concerning generalization for ML algorithms (SVMs, trees, Neural Networks) when the data is huge?


title: When to interpret multivariate tests when performing repeated-measures ANCOVA?

text: Yesterday, I ran a repeated-measures ANCOVA. The purpose was to determine the usability of two computer systems, call them 1 and 2. Each subject completed three (conceptually different) tasks, call them A B and C, on each system. (The covariate was whether the subject had prior experience with the current system, system 1). The dependent variable was time taken to complete each task. A GLM in SPSS prints out both the standard and multivariate results. The standard, within-subjects results indicated a main effect of system and an interaction of system*task. Now why would I need to use the multivariate results? They were significant as well, and I've heard multivariate tests have more power. But this seems to provide too much wiggle room for researchers-- if the within subjects test isn't significant, they can just look to see if the multivariate results are significant. Also, since there was an interaction system*task, I ran simple main effects on each task. Turned out, A and B were significant in one direction, C was significant in the other. This seems to suggest to me that doing a MANOVA is a bad idea here. But how would I even know that if I was only looking at the multivariate results? In short, I guess I just don't understand when to use MANOVA vs ANOVA if I have more than one DV.


title: Comparing two regression models

text: I have two regression models.   1. `lm(TEE ~ weight + gender)`   2. `lm(TEE ~ BMR)` How do I compare these two models and check which one predicts TEE best?


title: Why isn't there random effect?

text: I'm trying to use lme4::glmer to fit a mixed model such like that:               library(lme4)     set.seed(123)     df<-data.frame(id=sample(LETTERS[1:10], 50, T),                    y=rbinom(50, 1, 0.3),                    x1=rbinom(50, 1, 0.5),                    x2=as.integer(rnorm(50, 40, 5)))     df<-df[order(df$id),]          fitm<-glmer(y~x1+x2+(1|id), data=df, binomial)          coef(fitm)     $id       (Intercept)    x1       x2     A      -1.009 1.239 -0.01631     B      -1.009 1.239 -0.01631     C      -1.009 1.239 -0.01631     D      -1.009 1.239 -0.01631     E      -1.009 1.239 -0.01631     F      -1.009 1.239 -0.01631     G      -1.009 1.239 -0.01631     H      -1.009 1.239 -0.01631     I      -1.009 1.239 -0.01631     J      -1.009 1.239 -0.01631      I'm wondering why the effects are identical across ids. Why isn't there random effect as expected? Note, the glmer's out is exactly same as using `glm` here:               fitg<-glm(y~x1+x2, data=df, binomial)     coef(fitg)     (Intercept)          x1          x2         -1.00861     1.23882    -0.01631


title: What statistical model should I use to analyze the likelihood that a single event influenced longitudinal data

text: I am trying to find a formula, method, or model to use to analyze the likelihood that a specific event influenced some longitudinal data. I am having difficultly figuring out what to search for on Google. Here is an example scenario: Image you own a business that has an average of 100 walk-in customers every day. One day, you decide you want to increase the number of walk-in customers arriving at your store each day, so you pull a crazy stunt outside your store to get attention. Over the next week, you see on average 125 customers a day. Over the next few months, you again decide that you want to get some more business, and perhaps sustain it a bit longer, so you try some other random things to get more customers in your store. Unfortunately, you are not the best marketer, and some of your tactics have little or no effect, and others even have a negative impact. What methodology could I use to determine the probability that any one individual event positively or negatively impacted the number of walk-in customers? I am fully aware that correlation does not necessarily equal causation, but what methods could I use to determine the likely increase or decrease in your business's daily walk in client's following a specific event? I am not interested in analyzing whether or not there is a correlation between your attempts to increase the number of walk-in customers, but rather whether or not any one single event, independent of all others, was impactful. I realize that this example is rather contrived and simplistic, so I will also give you a brief description of the actual data that I am using: I am attempting to determine the impact that a particular marketing agency has on their client's website when they publish new content, perform social media campaigns, etc. For any one specific agency, they may have anywhere from 1 to 500 clients. Each client has websites ranging in size from 5 pages to well over 1 million. Over the course of the past 5 year, each agency has annotated all of their work for each client, including the type of work that was done, the number of webpages on a website that were influenced, the number of hours spent, etc. Using the above data, which I have assembled into a data warehouse (placed into a bunch of star/snowflake schemas), I need to determine how likely it was that any one piece of work (any one event in time) had an impact on the traffic hitting any/all pages influenced by a specific piece of work. I have created models for 40 different types of content that are found on a website that describes the typical traffic pattern a page with said content type might experience from launch date until present. Normalized relative to the appropriate model, I need to determine the highest and lowest number of increased or decreased visitors a specific page received as the result of a specific piece of work. While I have experience with basic data analysis (linear and multiple regression, correlation, etc), I am at a loss for how to approach solving this problem. Whereas in the past I have typically analyzed data with multiple measurements for a given axis (for example temperature vs thirst vs animal and determined the impact on thirst that increased temperate has across animals), I feel that above, I am attempting to analyze the impact of a single event at some point in time for a non-linear, but predictable (or at least model-able), longitudinal dataset. I am stumped :( Any help, tips, pointers, recommendations, or directions would be extremely helpful and I would be eternally grateful!!! Thank you


title: How to: repeated measures anova and Friedman's test with multiple imputations for missing data

text: I have a study where patients were measured pre-operatively (time_pre), post- operatively (time_post) and long-term follow-up (time_lt). I have 60 subjects. Data consists of: 5 scale variables (1-100) 14 ordinal variables (1-5) and some patient specific characteristics (none missing) like age, gender, co- morbidities (yes / no) I want to check if there is a significant difference between pre-operative scores and post-operative scores (long- and short-term) and between short and long-term follow-up. Therefore I used a repeated measures ANOVA on the scale variables and with post-hoc Bonferroni correction, and Friedman's test for the ordinal variables with post-hoc Bonferroni correction. Now I want to repeat these tests, but with the missing values filled in. Therefore I performed multiple imputations, but I can't figure out how I should perform a repeated measures ANOVA from this point. I use SPSS v22 statistical software.


title: Classifier or regression for binary system?

text: Lets say you want to create a random forest model that predicts whether a user will click on your search engine ad. Lets say the training data set marks the independent variable (i.e. whether the user clicked or not) as a 0 if it's a miss, and a 1 if it's a hit. Should this type of model be a classifier or a regression model? Or can it be either? **Background** I am trying this both ways, and I am finding that in the classifier version of this model, the model predicts a miss nearly every time. I'm not sure if there is a problem in the way I've defined the model, or if this is by design, because for a given data point, if you repeated that data point 100 times it might only be a hit 10% of the time, even for the likely-to-be-a-hit data points. In other words, there is almost no single point where the user is actually likely to click on the ad, so I think that maybe given that the probability is always less than 0.5, the model is simply predicting a miss every time


title: Two-stage linear regression

text: If I do regression in two stages: > Stage 1: $y\sim x_1 + 1$ > > Stage 2: resid_1st_stage $\sim x_2 + 1$ Will the `resid_2nd_stage` be orthogonal to $x_1$?


title: Can I mix different data transformations in the same model?

text: I assume it is OK to mix different data transformations in the same analysis. I've had to transform some variables to squared and some to cubed in order to meet normal distribution requirements. I'm assuming that it's OK to then use the transformed variables together for regression analysis?


title: Advanced statistics books recommendation

text: There are several threads on this site for book recommendations on introductory statistics and machine learning but I am looking for a text on advanced statistics including, in order of priority: _maximum likelihood, generalized linear models, principal component analysis, non-linear models_. I've tried Statistical Models by A.C. Davison but frankly I had to put it down after 2 chapters. The text is encyclopedic in its coverage and mathematical treats but, as a practitioner, I like to approach subjects by understanding the intuition first, and then delve into the mathematical background. These are some texts that I consider outstanding for their pedagogical value. I would like to find an equivalent for the more advanced subjects I mentioned.   * Statistics, D. Freedman, R. Pisani, R. Purves.   * Forecasting: Methods and Applications, R. Hyndman et al.   * Multiple Regression and Beyond, T. Z. Keith   * Applying Contemporary Statistical Techniques, Rand R. Wilcox


title: Prediction in Cox regression

text: I am doing a multivariate Cox regression, I have my significant independent variables and beta values. The model fits to my data very well. Now, I would like to use my model and predict the survival of a new observation. I am unclear how to do this with a Cox model. In a linear or logistic regression, it would be easy, just put the values of new observation in the regression and multiply them with betas and so I have the prediction of my outcome. How can I determine my baseline hazard? I need it in addition to computing the prediction. How is this done in a Cox model?


title: Methods of "Two-Way" Multidimensional Scaling (MDS)?

text: I just want to ask about the MDS. Below is the subcategory of the Multidimensional Scaling topic. **Multidimensional Scaling**   1. Metric and Non-Metric Models   2. Methods of "Two-Way" Multidimensional Scaling (MDS) I was tasked to report the second sub-topic of MDS which is the **Methods of "Two-Way" Multidimensional Scaling**. I don't have any idea about this topic actually, but when I try to surf through the internet, I found out one book which is all about **Multidimensional Scaling**. The problem is that I can't find the topic. But, I observed that there is this Section of Chapters entitled **MDS and Related Methods** consisting of chapters entitled _Procrustes Procedures_ , _Three-Way Procrustean Models_ , _Three-Way MDS Models_ , _Modeling Asymmetric Data_ , and _Methods Related to MDS_. I know this could not be the **Methods of "Two-Way" Multidimensional Scaling (MDS)**. I'm not sure. I just want a help on the description of the **Methods of "Two-Way" Multidimensional Scaling (MDS)** , or what does it talks about? or what are the possible sub-topics to tackle if any. Your help would greatly be appreciated.


title: Regarding sum of squares two ways how are they connected

text: $$ \sum(X^2) - \frac{(\sum X)^2}{n} = \sum(X^2) - m\sum X $$ This was nicely derived in Sum of squares two ways, how are they connected? But why is the second term called the "correction term for the mean"? This is lost on me. What is needing correction and why?


title: One-sided Fisher's exact test and its complement

text: If a one-sided Fisher's exact test testing the null hypothesis that A is not superior to B has a p-value of 0.98, is it true that the Fisher's exact test p-value for the null hypothesis that B is not superior to A is 1-0.98 = 0.02? The question arises from the paper at nejm.org/doi/full/10.1056/NEJMoa1100403 where the primary outcome had a one-sided p-value of 0.98, but when I calculated the P value with a one-sided FET, I got a p-value of 0.09. Therefore, wondering why the "other" one-sided FET wouldn't have a p-value of 1-0.09 = 0.91 (instead of 0.98). Would appreciate help in knowing where I'm going wrong. The Stata command used to generate the p-value of 0.09 was:                csi  2  7 2850 2854 , exact                           |   Exposed   Unexposed  |      Total     -----------------+------------------------+------------                Cases |         2           7  |          9             Noncases |      2850        2854  |       5704     -----------------+------------------------+------------                Total |      2852        2861  |       5713                      |                        |                 Risk |  .0007013    .0024467  |   .0015754                      |                        |                      |      Point estimate    |    [95% Conf. Interval]                      |------------------------+------------------------      Risk difference |        -.0017454       |   -.0037999    .0003091            Risk ratio |         .2866159       |    .0595926    1.378504       Prev. frac. ex. |         .7133841       |   -.3785042    .9404074       Prev. frac. pop |         .3561301       |                      +-------------------------------------------------                                       1-sided Fisher's exact P = 0.0904                                       2-sided Fisher's exact P = 0.1793


title: Bootstrap with replacement with small number of repetitions

text: In this youtube video about bootstrap resampling, the creator states that when the number of bootstrap processes is small, the distribution of the parameter being estimated can no longer be thought to have a normal distribution (see from the 8:50 mark). Since the parameter's distribution is not normal, the standard deviation can't be obtained via the standard formula: ![enter image description here](http://i.stack.imgur.com/63JGV.png) The author thus states that one needs to use the CDF to calculate where the 68% confidence interval values are and obtain from there an estimation of the standard deviation. My question is: where can I get a source for this workaround for obtaining the standard deviation when the number of times the bootstrap process is repeated is small? * * * I'm interested in this because I'm using bootstrap with replacement to estimate the uncertainty of a parameter and I can't possibly repeat it thousands of times (not even dozens actually) because it would be impossibly time consuming.


title: Early detection of peaks with neural network

text: I am using neural network technique (backward learning). As a output for example I am giving the 18 point ahead value and as input I gave the latest 5 point to train. (I tried the many combinations of input data 5, 10, 20, 30...). For example, the way I trained my data:               t, t+1, t+2, t+3, t+4... => t+22(4+18)     t+1, t+2, t+3, t+4, t+5... => t+23      Exponential inputs:               t, t+1, t+2, t+4, t+8...   => t+26(8+18)     t+1, t+2, t+3, t+8, t+9... => t+27      After I trained, I have done forward learning with the my trained values. I have observed that neural network won't able to catch the sudden peaks. Most of the time if I am going to predict 18 seconds ahead, it predicts the correct result 17 seconds later.   * In order to predict long term predictions for example 1000 node later on, what will be the best way to pick input values?   * Do you have any advice for me about how could I able to predict sudden peaks(that will happen t seconds later) with neural network? ![enter image description here](http://i.stack.imgur.com/6dq9W.gif)


title: Explicit mapping of input space to a high dimensional feature space

text: I am now studying on Kernel Methods, all the theory behind it to understand Support Vector Machines. Of course, I understood some very well but there is something I could not completely comprehend. Implicit/Explicit mapping of the input space to feature space is so important for kernel-related machine learning techniques. With explicit mapping, meaning that only using a mapping function $\phi$, it is possible that inputs can be projected to a infinite feature space or sometimes, it is hard to work on high-order polynomials with many coefficients, increasing dimensionality of feature space and number of free parameters So we have kernel trick, it defines a similarity measure between two vectors in an input space from the dot product $k(x_i, x_j)=<\phi(x_i), \phi(x_j)>$. So that one can work on implicitly on a high dimension without doing lots of computations. That's all ok for now. The point I do not understand is about explicit mapping, i.e. only using $\phi$ function. I was just looking for an example online, I run into with that one below. ![enter image description here](http://i.stack.imgur.com/SVYgd.jpg) For example, here we have circular decision boundary given in $\mathbf{R}^{2}:x_{1}^{2}+x_{2}^{2}=1$.Thus, data is replotted in a $\mathbf{R}^{3}$ feature space using quadratic map $\phi(x_{1}, x_{2})=(x_{1}^{2}, x_{2}^{2}, \sqrt{2}x_{1}x_{2})$ and is then linearly separable, thus we can work with linear methods on a transformed data. This approach is called preprocess and is mostly used from the beginning of Machine Learning. By selecting $ x_{1}^{2}$ , $ x_{2}^{2}$ , $ \sqrt{2}x_{1}x_{2}$ as features, the new feature space can be constructed and for separation, all we need is a linear hyperplane. > But, what about the case if we need to map an input space which has 70 x 70 > pixels as patterns to a feature space. In this example, if we choose as > nonlinearity all 5th order monomials and since our dimension vectors have > 4900 rows, the dimension of the feature space can be computed as below. > > $\left(\begin{array}{c} {5+4900-1} \\\ {5} \end{array}\right)\approx > 2.3\times 10^{17}$ > > To make such a mapping is almost impossible. I cannot understand this example in quotation. First of all, how is my data with 70X70 pixels and how is it possible to have 4900 rows? I cannot imagine the data set in my mind. What does he mean with "nonlinearity all 5th order monomials"? I am not good at image classification so this might be the reason. I know it will be sound a bit silly but I will be much obliged if you care to explain it in details. Thanks in advance.


title: SVM Laplacian kernel parameters with R

text: In this paper, the authors write the SVM Laplacian kernel as follows: ![enter image description here](http://i.stack.imgur.com/pOJWO.png) I would like to try the same and tune the 3 parameters ("p",a,b). I am using the R kernlab ksvm function (manual here, page 54-57), but I do not know how to match the parameters to the one mentioned in the formula above. Could anyone advise which ksvm parameters match ("p",a,b) above ? Or perhaps another R library for SVM?


title: What sense does it make to compare p-values to each other?

text: I have two populations, each containing 1000 samples. For each sample I have two properties A & B. I have used a t-test separately for A & B: both were found significantly different within the two groups; A with P=0.008 and B with P=0.002. Is it okay to claim that the B property is better discriminated (more significant) then the A property? Or is it that a t-test is just a yes or no (significant or not significant) measure? **updated** : The measurements for two different properties 1k each, men and women, (A) first-year grade point average and (B) SAT scores. **updated2** : according to the comments here and to what I have read in p-value. I think that the answer should be: drop the meaningless p-value and reports your effect size. any Idea ?


title: Is it ok to transform a logarithm variable to z score

text: I have a variable that has 57 kurtosis, so I decided to transform it to log. However, I have multicolleanirity problem due to interacting this variable and others with another variable so I am using z scores to reduce VIF value. So ... > > Is it ok to transform the log variable to z score thank you so much


title: Estimating the probability of next failure?

text: I have some data about how some components failed. They can either be replaced in whole or have some of their sub-components repaired. I would like to know if there is a way to fit a failure model to the underlying failure distribution and computing reliability of my components. I read about Survival Analysis and was very interested in it (because it also takes into account the components that have "survived" post observation period), but unfortunately, it assumes complete component failures (whereas my situation deals with component-repairs as well). Is there any related area that deals with the repair situation? The data I have is the following (output of `dput` in `R`):               c(43044461L, 48043852L, 39290306L, 50335198L, 47604418L, 9875530L,      1738141L, 46301014L, 6477432L, 7372344L, 49950268L, 29007573L,      46118058L, 12292465L, 238234L, 9277711L, 34570342L, 37371531L,      28471868L, 11157523L, 50931534L, 26694512L, 29961599L, 2457994L,      33987131L, 20175414L, 1577647L, 769861L, 1249158L, 32224970L,      9409192L, 34033673L, 32297668L, 27564068L, 2730862L, 19241848L,      18254820L, 39185109L, 33662797L, 12003705L, 39384318L, 26484711L,      15432679L, 32281540L, 34839587L, 385860L, 9973362L, 21839720L,      37642093L, 33886689L, 17704281L, 33337960L, 13913748L, 32304197L,      33636487L, 22603695L, 25718240L, 32181213L, 32264085L, 5881735L,      35120576L, 5527473L, 23051956L, 23054083L, 22961394L, 22962199L,      4632371L, 2757599L, 9663000L, 29225898L, 29233148L, 29237159L,      29435753L, 1107666L, 26138387L, 26163543L, 26215242L, 7561251L,      26218513L, 26221383L, 30476747L, 26222109L, 27602048L, 27607598L,      27610958L, 27684885L, 1782502L, 25654474L, 19616668L, 10155706L,      1736017L, 27025910L, 10117992L, 20506415L, 11541948L, 27586115L,      26692258L, 16235808L, 14966802L, 15974788L, 2014068L, 7215767L,      8022835L, 22321409L, 22524449L, 9940245L, 12431361L, 6723321L,      1337135L, 14011887L, 26516990L, 19894509L, 23580559L, 25837938L,      25421604L, 25260751L, 1014468L, 17600421L, 18969565L, 19501830L,      11816652L, 11832654L, 18962013L, 18967493L, 19044885L, 10660220L,      7152965L, 25676536L, 182572L, 14811759L, 18464712L, 12164906L,      19920052L, 12495744L, 16548293L, 16625312L, 9665000L, 284178L,      284213L, 284143L, 284178L, 284178L, 284178L, 2992895L, 9475072L,      4444603L, 10436184L, 3108530L, 8261921L, 3493839L, 8603179L,      563679L, 7788458L, 12934577L, 16212035L, 14764962L, 3274975L,      18208633L, 2661185L, 12233593L, 19688931L, 3176516L, 10347964L,      3176445L, 7326030L, 816261L, 17698452L, 8280151L, 9832881L, 18059503L,      16362051L, 18559410L, 3016907L, 18821969L, 250538L, 18937854L,      20460554L, 11466655L, 19841035L, 5666426L, 11466795L, 9351749L,      14202233L, 8395799L, 2469799L, 11120915L, 2007822L, 2251827L,      10425100L, 12864547L, 17400937L, 12224074L, 14425083L, 7790951L,      7790951L, 18092421L, 249048L, 14389608L, 14389643L, 14389643L,      14303313L, 14303313L, 184829L, 184797L, 184972L, 13871348L, 13871243L,      13871383L, 13871278L, 99619L, 99514L, 98467L, 98467L, 98429L,      98467L, 98467L, 98432L, 98432L, 98432L, 98432L, 98432L, 11372596L,      9926436L, 204260L, 2492418L, 2492418L, 2960708L, 9660687L, 8180564L,      53600L, 11428758L, 9942851L, 2648374L, 1899270L, 1812870L, 1812870L,      1812870L, 1812870L, 1812870L, 1812870L, 5811703L, 3334565L, 11639886L,      5776214L, 5386427L, 2057170L, 853094L, 4856466L, 7488042L, 10322886L,      71544L, 2208165L, 2035661L, 2121240L, 2121744L, 2122268L, 4393081L,      2481686L, 4767429L, 2696190L, 7631815L, 3464979L, 2737688L, 11709771L,      9311182L, 6952803L, 210747L, 8219827L, 508306L, 4713143L, 6456962L,      5299429L, 8969430L, 2290297L, 117516L, 2980892L, 2048608L, 7152852L,      396062L, 396062L, 3459981L, 3691834L, 3834119L, 6950281L, 1431240L,      15124240L, 1292511L, 14712416L, 14629719L, 15120090L, 3284021L,      1485496L, 2420495L)


title: Decrease in p-value when using student-t instead of normal distribution for GARCH?

text: I am trying to model the return and volatility of different financial data such bmw(evir package in R),SP500(MASS package in R),gold and oil. Their statistics indicate departure from normality and the Jarque Bera test reject the hypothesis of normality. I tried using the student t instead and got a lower AIC which I suppose is to be expected. However, the p-value for standardised squared residuals significantly decreased sometimes even below 0.05. In addition when I tried an out-sample test the normal model out- performed the student-t model in all cases. The MSE for both was around 2 for both models which I think is very high as I saw a study getting MSE around 0.08. I'm not sure what is wrong here as I believe the student-t should be the one out-performing and I'm shocked with such a high MSE. I used ARMA(p,q)/GARCH(1,1) where p=0,1 and q=0,1 but none gave me the results that I was hoping for.


title: Statistical test for a significant change in time series (sales) trend after policy change

text: My apologies if this has been addressed before, however my experience with statistics is a bit limited, particularly when it comes to time series analyses. My current question relates to the analysis of a trend in monthly sales before and after an event for a single product. Suppose there is a new product with increasing monthly sales, however the company decides to initiate a promotion 12 months after it launches to further increase sales. A year later, the company wants to examine if the sales trend after the promotion is significantly different from the trend prior to the promotion (ideally the unit sales would be increasing at a faster rate after initiation of the promotion). So essentially, I want to determine if I had forecasted 12 months into the future just prior to the promotion, are the sales significantly different from what they actually turned out to be. What statistical test is the most appropriate to use to determine if the monthly sales trend is significantly different before vs. after the promotion? I only have one data point for each time point (monthly sales) before and after, and since the sales are increasing I have come across problems due to the "stationarity" requirement for many tests, as well as the fact that I would only have 12 historical data points on which to base the initial forecast. Any help would be greatly appreciated. Thanks.


title: How to log transform data with a large number of zeros

text: I have a dataset comprised of continuous values that have about 30-50% zeros and a large range (10^3 - 10^10). I believe these zeros are not a result of missing data and are the result of the sensitivity of the machine taking the measurements. I would like to log10 transform this data so I can look at the distribution, but I'm not sure how to handle the zeros I've done a lot of searching and found the following   1. Add a small constant to the data like 0.5 and then log transform    2. something called a boxcox transformation I looked up boxcox transformation and I only found it in regards to making a regression model. I just want to visualize the distribution and see how it is distributed. Currently when I plot a historgram of data it looks like this ![enter image description here](http://i.stack.imgur.com/OPraf.jpg) When I add a small constant 0.5 and log10 transform it looks like this ![enter image description here](http://i.stack.imgur.com/8Nzu0.jpg) Is there a better way to visualize the distribution of this data? Im just trying to get a handle on what the data looks like in order to figure out what kind of tests are appropriate for it.


title: How to do multivariate regression in R?

text: I need to preform a multivariate normal regression in R. The question is: > Let $Y_1$, $Y_2$, and $Y_3$ follows multivariate normal distribution. What > is > >   1. the conditional of $Y_3$ given $Y_1$ and $Y_2$ >   2. the conditional of $Y_2$ given $Y_1$ > > From these two, derive: > >   3. the joint distribution of $Y_3$ and $Y_2$ given $Y_1$. > > Now suppose you have a sample of size $n$ from the multivariate normal distribution. Do the two regressions in (1) and (2). **How can I combine them to get (3), the regression of $Y_3$ and $Y_2$ on $Y_1$?**               library(mvtnorm)     mu  <- c(1,2,3)     Sig <- matrix(c(4,2,1,2,4,-1,1,-1,4), nrow=3, ncol=3)      Y   <- rmvnorm(20, mean=mu, sigma=Sig) #generate multivariate normal distribution     y3  <- lm(Y[,3]~Y[,1] + Y[,2])      y2  <- lm(Y[,2]~Y[,1])


title: How to convert log likelihoods into scores in Naive Bayes?

text: I am currently implementing a text classification program with Naive Bayes. I produce two multinominal models in my training function: p(w|nonSPAM) and p(w|SPAM)) as well as a prior probability P(S). In my testing function I go through each test document, and for each test document, I go through all the terms and compute logP(nonSPAM|D) and logP(SPAM|D). Then I make a classification decision by comparing these two quantities (SPAM = 0 or nonSPAM =1). My problem is: I want to return a score (e.g. 0.52121) rather then a strict probability (0 or 1), so I can use different thresholds in my program. Is it possible to calculate a score only by using logP(nonSPAM|D) and logP(SPAM|D) and prior probability P(S)? I asked a similar question here but my current question is more related with statistics.


title: Approximation to Lognormal Distribution

text: I'm looking for an approximation to the curve of a lognormal distribution, for use in non-linear regression against a dataset. As an alternative, I'm interested in an approximation to the CDF thereof. I have several goals:   1. Determine how closely a sample dataset drawn from a process of unknown distribution matches a lognormal distribution.   2. Given a set of samples drawn from a process determined or known to be lognormally distributed, but with unknown characteristics, determine the probability of a observing a future sample with a given value that may lie outside the range of values observed thus far.   3. Given a set of samples as above, and a new sample that may lie outside the range of values in the samples already seen, determine the likelihood of having observed that sample.   4. While computing the above, efficiency and simplicity of the implementation is quite important. A fast approximation with well-understood properties such as error bounds is better than an accurate algorithm that is difficult to implement or computationally intensive. For these purposes, I think it is best if the approximation has a closed form so that a nonlinear regression can be done with lower computational overhead. Ideally if it is an approximation to the lognormal distribution itself, then it would be nice if it has a simple integral as well so that the CDF can be approximated too. It's possible that I'm trying to go about this the wrong way. Here's a more specific question to help figure that out: suppose I gather 1000 samples from my process. The values of the samples mostly range from 1 to 10, with occasional samples up to 20 or so. I know (based on experience) that long-term it is quite possible to see samples in the range of 10 times higher than that, but I haven't observed any from _this_ process yet. How can I determine the probability of the 1001st sample having a value greater than 100 or another arbitrary number? If the 1001st sample's value is 180, how can I determine how likely that was, based on the first 1000 samples?


title: Anomaly/outlier detection using fuzzy clustering

text: I understand that fuzzy clustering using FCM produces a membership matrix for the set of data points we feed to it. What characteristics will an anomalous cluster produced during this method have? (Considering I only have unlabelled data)


title: Explanation of "R Square" in terms of correlation?

text: In linear regression, I seemed to see (but I forgot where) that if the square sum of residuals divided by square sum of total variation equals the square of some correlation between some random variables, and hence it explains why the ratio is called "R Square". So I was wondering what the correlation is between, and why the ratio is the square of the correlation? What assumptions on the linear regression model are needed for the explanation of "R Square" in terms of correlation? Thanks!


title: Random sampling (Real data) is so important, why?

text: How it is understandable for computer (or software) that a sampling method is Random? In fact, it is possible to data sampling by random method in N cases but one of the N cases been like as non-random! That means sampled data randomly is like non-randomly. Now, in this case, how computer find out that selecting samples distribution is based on random selecting? and why random sampling is so important? * * * ### Edit Ok. thanks for answer but also i have a question: Assume your sampled data located at almost one case e.g. age in heigh of population. so, you have a pattern in collected data while i did random sampling (this is one of cese of random sampling). now, what is important in random sampling than non-random that sampled data is accepted to following analysis? why random sampling is important to non-random?


title: How to find out if a set of daily measurements  are random or not?

text: There is a set of daily measurements. Time and measured values are both discrete. I want to find out whether measured values depend on the day the measurement was taken, or whether measurements are completely random. In other words, I want to find out if it is possible to predict measured values of a certain day or not.   * What subjects of statistics should I study to be able to solve this problem?  Please give me some keywords or direction. **Edit**   Some more information. Imagine the following situation. A machine chooses each day a letter (simply a byte) and displays it on a screen. The process which is used to choose daily letter is unknown, but it is clearly algorithmic (not measure the wind speed or count people in the room or similar). Someone collected part of daily letter over a period of time and now wants to understand if it is possible to produce the next (or any) daily letter. Some methods possibly employed by the machine are considered "hard" (or random). For example use a secret key to encrypt the current date (predicting the next letter will be equivalent in most cases to braking the encryption). Other are considered "easy", for example xor of all bytes in date's representation. If the method is random, then there is no hope to produce the next daily letter.


title: Summarizing k-fold cross validation results

text: I have a project for which I had to do a 5-fold cross-validation. The dataset comprises of 5 days, so we take one day as our training set and apply on the remaining 4, rinse and repeat till we get full coverage. Pretty standard stuff. Now, here's my question, once I finished the cross-validation, I am left with True Positive and False Positive rate for each one my runs. Is there a standard way to show this data visually that would summarize the results ? In particular I want to show that the technique is fairly stable and about a day's worth of data is enough to train it successfully.


title: Getting a zero p-value with Wilcoxon test in Matlab

text: I'm creating a classification model for detection of stress. I'm trying to select the best features from 56 total. I have a matrix of `1937x56` with the non stress data, and a matrix of `1763x56` for stress data. `f_n -> non-stress matrix` and `f_s -> the stress one` I'm using `Matlab` and here's my code.               rank = zeros(1,size(f_n,2));          for l = 1:size(f_n,2)         rank(1,l) = ranksum(f_n(:,l), f_s(:,l));     end      the think is I'm getting some `p-values` equal to zero. Here's the output of `rank`               3.23970807719953e-59    1.19012268294050e-32    1.24342293425171e-123   2.98338574026435e-227   1.94090447559927e-227   7.54664151060207e-250   8.99703327997911e-292   1.38826281436916e-298   4.78404773111577e-308   2.15213684390080e-37    0.000196579371678330    4.25370029962320e-27    4.90973329410252e-166   9.84016939125288e-200   4.84716291687399e-277   3.52745590178140e-287   2.81025758599302e-297   7.01405396494596e-289   3.30447666181206e-165   1.73163324556308e-32    1.06353645618514e-116   6.74904729244135e-130   3.77342641391818e-173   4.64192483365192e-239   1.60123198146693e-250   1.16096192265296e-259   1.31621596229152e-271   0.00115260954904737 5.11414519751266e-253   0   0   0   0   0   0   0   5.53879714415676e-137   5.61283198666988e-67    2.53913294535256e-177   0   0   0   0   0   0   1.70864676934011e-91    4.50477630440800e-73    1.34942550217303e-165   4.91843135947255e-180   8.45438406534899e-209   1.02062206403179e-213   1.26764781876251e-214   6.91325407956207e-247   2.35709504234097e-276   0   1.11140566121933e-36      I don't think it is supposed to be really zero. Any help?


title: Generating values from a multivariate Gaussian distribution

text: I am currently trying to simulate values of a $N$-dimensional random variable $X$ that has a multivariate normal distribution with mean vector $\mu = (\mu_1,...,\mu_N)^T$ and covariance matrix $S$. I am hoping to use a procedure similar to the inverse CDF method, meaning that I want to first generate a $N$-dimensional uniform random variable $U$ and then plug that into the inverse CDF of this distribution, so to generate value $X$. I am having issues because the procedure is not well documented and there are slight differences between the mvnrnd function in MATLAB and a description that I found on Wikipedia. In my case, I am also choosing the parameters of the distribution randomly. In particular, I generate each of the means, $\mu_i$, from a uniform distribution $U(20,40)$. I then build the covariance matrix $S$ using the following procedure:   1. Create a lower triangular matrix $L$ where $L(i,i) = 1$ for $i=1..N$ and $L(i,j) = U(-1,1)$ for $i < j$   2. Let $S = LL^T$ where $L^T$ denotes the transpose of $L$. This procedure allows me to ensure that $S$ is symmetric and positive definite. It also provides a lower triangular matrix $L$ so that $S = LL^T$, which I believe is required to generate values from the distribution. Using the guidelines on Wikipedia, I should be able to generate values of $X$ using a $N$-dimensional uniform as follows:   * $X = \mu + L * \Phi^{-1}(U)$ According to the MATLAB function however, this is typically done as:   * $X = \mu + L^T * \Phi^{-1}(U)$ Where $\Phi^{-1}$ is the inverse CDF of a $N$-dimensional, separable, normal distribution, and the only difference between both methods is simply whether to use $L$ or $L^T$. Is MATLAB or Wikipedia the way to go? Or are both wrong?


title: Mann-Whitney test with ties

text: I have two data sets with PDFs roughly like this: $$ p(x) = \left\\{ \begin{array}{lr} .75 & x = 0\\\ \text{Lomax}(x) & x > 0 \end{array} \right. $$ i.e. it's continuous, except there is a huge number of points at $x=0$. I wish to tell if the values in one set are lower than the other. However there are a lot of ties in my data, which makes it kind of hard. It was suggested to me that I do something like: bootstrap a bunch and if a 95% confidence interval has x's being less than y's (as determined by Mann-Whitney U), we declare $X<Y$. It's not really clear from the definition of U what I do about ties though. `wilcox.test` in R seems to break them evenly, e.g. when comparing $(0,0)$ to $(0,0)$ it says $U=2$. (This behavior seems to contradict its documentation though which says it returns "the number of all pairs (x[i], y[j]) for which y[j] is not greater than x[i].") Is breaking them evenly the right thing to do here? I'm worried since I have such a huge fraction of ties that my intuition about a "95% confidence interval" doesn't actually correspond to 95% of my samples.


title: Time series data distribution forecast?

text: While having chronically data of population growth (registered users of a site), I want to compute a function that approximates future growth, based on past data. Also, what we ll be the distribution of that function? What is the distribution of interarrivals between consecutive registrations.


title: More complicated classifier vs data preprocessing

text: For example if we have two options to use non-linear classifier like SVM with kernel or use linear classifier like linear SVM with data preprocessing like some non-linear dimensionality reduction which one is better? In other words why we should use a complicated classifier if we can do some data preprocessing (or "projection" if I am not mistaken it's called manifold learning) to make data linearly separable?


title: About comparing standard scores

text: I have time-series data, from which I would like to detect anomalies. I have Z-scored the time-series, and there is two points which clearly deviate from the baseline. Their Z-scores are 3 and 6, respectively. So, I can say they differ significantly from the baseline. Can I also say that the two anomalies themselves differ significantly, given their different Z-scores? Needless to say I'm new to statistics...


title: How do I transform my data so that it has mean zero and standard deviation one?

text: Given a series of values,               3.00,5.00,7.00,4.00,7.00,5.00,22.00,4.00,6.00,7.00,9.00,6.00,4.00         I need to 'rescale' the data so they have new values with a mean of 0 and a standard deviation of 1. Having followed a number of examples for the equation, I get the following,               Mean    6.85     StdDev  4.85      The standard deviation I worked out the long way via taking an original value, subtracting the Mean, squaring the result, summing all the squares, dividing that sum by the count -1 before finally getting the square root of that.               (3.00 - 6.85) = -3.85    >    -3.85² = 14.79     (5.00 - 6.85) = -1.85    >    -1.85² =  3.41     ...     (4.00 - 6.85) = -2.85    >    -2.85² =  8.10          Sum                                   281.69     Variance (using n -1)                  23.47     Standard Deviation                      4.85      I'm sure it's correct as I checked it here http://www.mathsisfun.com/data/standard-deviation-calculator.html Using the mean and standard deviation I get new values of,               -0.79,-0.38,0.03,-0.59,0.03,-0.38,3.13,-0.59,-0.17,0.03,0.44,-0.17,-0.59      These values were calculated in the following manner               ((3.00 - 6.85)/4.85) = -0.79     ((5.00 - 6.85)/4.85) = -0.38     ....     ((4.00 - 6.85)/4.85) = -0.59      So, to my question. It's the standard deviation that is the confusing part. I want my new values to have a Mean of 0, and a standard deviation of 1. When I work out the SD for my original values, I get an SD of 4.85. What do I need to do to get an SD of 1 ? Thank you for taking the time to read the question. John.


title: Finding the highest probability balls

text: I'd appreciate any insights or references to research regarding the following: Suppose you have a discrete metric space with a probability distribution on it. also, suppose that I'm given a number $k$. My problem is finding the highest probability ball. A **ball** is the set of all points whose distance from a given point is less than or equal to $k$, and the "probability of the ball" is just the sum of the probabilities of this set of points. The continuous case (less interesting in my case, but potentially helpful) would be to have a probability distribution on $R^n$, an finding the n-ball with the highest integral when integrating the density function on the $n$-ball. This problem is an abstraction of a practical problem I try to solve. the solution of special case of $k=0$ (e.g. in which case a "ball" is just a single point) is just the mode of the distribution. The general case arises for example in the following prediction problem: suppose you try to "predict" $n$ independent binary events. So every prediction is a binary $n$-tuple (with a probability distribution imposed on that space, based, for example, on assumptions, previous data, etc.). Now, suppose you are allowed to have a few mistakes in your prediction, and you try to find a prediction that maximizes the probability of "success" (that is, not making too many mistakes). The solution is exactly an n dimensional "ball" in the n-tuple space with the hamming metric. Any idea, insight, or reference to relevant research would be much appreciated.


title: How to estimate 5 percent breaking point  of  a sample of ropes?

text: The breaking strengths in pounds of five specimens of rope were 660, 460, 540, 580, and 550. How can I estimate the point at which only 5 percent of such specimens would be expected to break?


title: Comparing models on the test data

text: My understanding is that if you (1) have a sufficiently large test dataset, and (2) your models have the same likelihood (noise assumption), then you should compare/select the model likelihood (or log-likelihood) on the test data, unadjusted by AIC/BIC/DIC/etc. Is there any justification for adjusting your likelihood for the test (or cross-validation) dataset for model selection? What is the best practice when such test dataset is available?


Is it possible to correct for violating the assumption of independence for nonparametric tests? I have a categorical independent variable and a categorical and binary dependent variable, and each subject was exposed to multiple levels of treatment. I have done some preliminary analysis with binary logistic regression and post-hoc chi-square tests, but I think these are inappropriate unless I can correct for my non-independent data. Info on this has been surprisingly difficult to find, and help is much appreciated.


title: Binomial Probability Question

text: Given two basketball players. John made 38/50 free throws. Mike made 80/100 free throws. What is probability that Mike is better at free throws than John?


title: Are sample means for quantiles of sorted data unbiased estimators of the true means?

text: In studying income inequality, it is very common to look at sample means for deciles or quintiles of the sample, and to assume that the sample means are good estimators of the true means. In this setting, the "deciles" and "quintiles" normally refer, not to the break points, but to the sets of observations divided by the break points. Suppose that the income values are observed with error, and that the error, or more probably, the percentage error, is independent of the true value.   * Is the mean of the sample quantile, e.g., the top decile, an unbiased estimator of the population mean? I know that with some leptokurtic distributions (e.g., the Pareto) the sample mean understates the population mean. My question refers, not to this, but to any bias that might be induced by the sorting process, because one sorts on the observed values including error, rather than on the true values.   * My intuition is that sample mean of the highest decile/quintile would be biased upwards, because positive errors would be sorted up, and conversely for the lowest. For instance, it seems to me that if incomes were non-negative but observed with a normal error, a large sample would contain some negative values, and a fine enough quantile would gather these negative values together into a lowest group with a negative mean, demonstrating bias since the true mean must be positive. Is this true?    * If the means are biased, is there any good way of correcting for this bias? Am I correct in thinking that the error's independence of the true value does not carry over to independence of the error from the _observed_ value, which includes the error? If so, is there an easy way to at least characterize, and ideally correct, that dependence?   * A commonly used index of inequality is the ration of mean of the incomes in the highest quintile or decile to the mean of the incomes in the lowest. If these means are biased, and the bias is corrected, will the resulting ratios be unbiased estimators of the true ratios?


title: Maximum value of coefficient of variation for bounded data set

text: In the discussion following a recent question about whether the standard deviation can exceed the mean, one question was raised briefly but never fully answered. So I am asking it here. Consider a set of $n$ nonnegative numbers $x_i$ where $0 \leq x_i \leq c$ for $1 \leq i \leq n$. It is not required that the $x_i$ be distinct, that is, the set could be a multiset. The mean and variance of the set are defined as $$\bar{x} = \frac{1}{n}\sum_{i=1}^n x_i, ~~ \sigma_x^2 = \frac{1}{n}\sum_{i=1}^n (x_i - \bar{x})^2 = \left(\frac{1}{n}\sum_{i=1}^n x_i^2\right) - \bar{x}^2$$ and the standard deviation is $\sigma_x$. Note that the set of numbers is _not_ a sample from a population and we are not estimating a population mean or a population variance. The question then is: > What is the maximum value of $\dfrac{\sigma_x}{\bar{x}}$, the coefficient of > variation, over all choices of the $x_i$'s in the interval $[0,c]$? The maximum value that I can find for $\frac{\sigma_x}{\bar{x}}$ is $\sqrt{n-1}$ which is achieved when $n-1$ of the $x_i$ have value $0$ and the remaining (outlier) $x_i$ has value $c$, giving $$\bar{x} = \frac{c}{n},~~ \frac{1}{n}\sum x_i^2 = \frac{c^2}{n} \Rightarrow \sigma_x = \sqrt{\frac{c^2}{n} - \frac{c^2}{n^2}} = \frac{c}{n}\sqrt{n-1}.$$ But this does not depend on $c$ at all, and I am wondering if larger values, possibly dependent on both $n$ and $c$, can be achieved. Any ideas? I am sure that this question has been studied in the statistical literature before, and so references, if not the actual results, would be much appreciated.


title: How to update a dynamic language model dataset?

text: I'm a statistics novice and I need help with a natural language problem. I'm writing a word-prediction algorithm for a mobile app. I'm using a unigram language model of word/count pairs where _count_ is the number of times that _word_ appears in a corpus. The algorithm is pre-loaded with a set of words and probabilities, _but it must also update the probabilities and learn new words as the user types_. Here is a sample of the pre-loaded table:               word                count     ----------        -------     you               1222421     I                 1052546     the                823661     to                 770161     a                  563578     and                480214     ...     spangles                5     moustached              5      The prediction algorithm is simple enough (rank predictions by descending _count_ ). My question is **how to update this table so the algorithm learns _quickly_ as the user types**. If I simply add new words with a count of 1, and increment count by 1 with each usage, the algorithm will learn too slowly. For example, suppose the user types "to" more often than "the". It will take over 50,000 usages before "to" outranks "the". That's slow learning! Similarly, suppose the user's name is "andy" and so he types that word more often than anything else. He'll have to type "andy" nearly half a million times before it outranks "and"! **How do I compute new probabilities so that so that it doesn't take hundreds of thousands of word usages to make a difference?** P.S. Computational efficiency is very important since this is for a mobile app. Sorry if I'm asking for too much! Thanks in advance!


title: How do you interpret a changepoint graph?

text: What is the correct way of interpreting a changepoint graph? Here is what R gives for changepoint analysis of the mean of a data set:               ----------     Changepoint type      : Change in mean      Method of analysis    : AMOC      Assumed Distribution  : Normal      Type of penalty       : SIC with value, 4.770685      Maximum no. of cpts   : 1      Changepoint Locations : 15 118       So changes occurred at t=15 and t = 118?


title: Null distribution on synthetic data

text: I need a sort of review of all (or most of) available methods to create null distributions (as reference) to use to compare a result. For ex: if I want to validate a network I can create a set of 1000 random networks starting from the initial dataset, by shuffling 1000 times my data. Then I obtain a null distribution of the 1000 nets toward which I can compare my own net (bootstrap) to score the significance (if my net is in the tails of the null distributions). Since this is not the best solution according to the data I start with, I need to know all the available methods used to create null distributions to use as a reference. Any suggestion about papers, links etc. useful to me to choose the best method? PS: I have to create the null distribution starting from a synthetic dataset (a dataset totally created in vitro) so I need to know how to build it first. Thanks in advance E.


title: Is there a term for min + (max - min) / 2?

text: I'm looking for a function to do that in R. I know how to write that function. Just don't want to reinvent something.


title: Framing the negative binomial distribution for DNA sequencing

text: The negative binomial distribution has become a popular model for count data (specifically the expected number of sequencing reads within a given region of the genome from a given experiment) in bioinformatics. Explanations vary:   * Some explain it as something that works like the Poisson distribution but has an additional parameter, allowing more freedom to model the true distribution, with a variance not necessarily equal to the mean   * Some explain it as a weighted mixture of Poisson distributions (with a gamma mixing distribution on the Poisson parameter) Is there a way to square these rationales with the traditional definition of a negative binomial distribution as modeling the number of successes of Bernoulli trials before seeing a certain number of failures? Or should I just think of it as a happy coincidence that a weighted mixture of Poisson distributions with a gamma mixing distribution has the same probability mass function as the negative binomial?


title: Relationship between Intimal thickness (as assessed by score) and age or hormonal status

text: Is there a statistical program to calculate power for data generated by assigning a score (not necessarily linear ie score of 2 may not be equal to twice a score of 1)


title: Why autocovariances could fully characterise a time series?

text: I read from textbook that ' _Autocovariance can fully charaterise the time series_ ' joint distribution', I do not fully understand the connection between covariance and joint distribution here. Please explain to me, someone?


title: Formula for comparing the youth of different groups?

text: For purely personal reasons, I recently tried to write a Ruby script that pulled information about various school rosters, counting the numbers of freshmen...seniors on the team and using that data to spit out some information about the teams relative youth based on the amount of eligibility left. As far as I am concerned, actual age is not an issue here. The forumla I used was something simple:               eligibility = (FR * 4) + (SO * 3) + (JR * 2) + SR     average     = eligibility / (FR + SO + JR + SR)      Anyway, this is obviously a really simple formula. It doesn't put any weight on any one factor and the total eligibility number really doesn't tell me anything. I was wondering if anyone had any better ideas when it came to realizing a formula for this issue. I'm not sure if this is the forum for this kind of question, but I became interested in the solution the other day and thought this might be a good place to look for an answer.


title: How to pull data formatted for Matlab into R?

text: I just ran a vendor-supplied program that converts data from their proprietary binary format to what they call a "Matlab format". There are two resulting files from the conversion, a "filename.dat" and a "filename.m". The "filename.m" file seems to have column heading information in it. The second "filename.dat" file is >13Gigs in size and I can't load it into anything at the moment to inspect it. How do I go about loading this "Matlab format" data into R?


title: Constrained linear regression through a specified point

text: I have a point (x,y) that I need a linear regressor to pass through given a data set (X,Y). How do I implement this in R?


title: Hypothesis Testing applied to real life

text: I was wondering if it is possible to apply the method of hypothesis testing to real life. For example if someone can use it for decision making. I have always used this method for homework problems but maybe we can use this method as an aid in decision making. Therefore we could somehow know for example the probability of rejecting wrongly an alternative decision. What do you think on that? And it would be nice if someone can give an example if he thinks that this can be done.


title: Assessing DNA sequencing quality

text: I work with a lot of bar charts. In particular, these bar charts are of basecalls along segments on the human genome. Each point along the x-axis is one of the four nitrogenous bases(A,C,T,G) that compose DNA and the y-axis is essentially how many times a base was able to be "called" (or recognized by a sequencer machine, so as to sequence the genome, which is simply determining the identity of each base along the genome).   Many of these bar charts display roughly linear dropoffs (when the machines aren't able to get sufficient read depth) that fall to 0 or (almost-0) from plateau-like regions. Is there a straightforward algorithm for assigning these plots a score that reflects this tendency? Is stdev a good place to start? I'm a programmer but not much of a statistician.


title: Generating random vectors with constraints

text: I need to create random vectors of real numbers a_i satisfying the following constraints:               abs(a_i) < c_i;           sum(a_i)< A;        # sum of elements smaller than A     sum(b_i * a_i) < B; # weighted sum is smaller than B      aT*A*a < D          # quadratic multiplication with A smaller than D          where c_i, b_i, A, B, D are constants.      What would be the typical algorithm to generate efficiently this kind of vector?


title: Calculating the probability that demographic trends happened by chance

text: I ran an ad campaign on Facebook and I have impressions (ad-views), clicks, and sign-ups broken down by gender and age. Certain demographics have a higher rate of clicks and sign-ups when compared to their respective percent of impressions. What is the best way to determine if these differences represent significant results instead of random noise? I am trying to see how age and gender predict sign-up.


title: Multiple imputation for clustered data

text: I have a few questions regarding multiple imputation for nested data. Context: I have repeated measures (4 times) from a survey and these are clustered in workplaces (205 workplaces). There are about 180 items on this survey. q1. Is it possible to take both the repeated measures and the workplace clustering into consideration or do i have to decide for one of the two? q2. If i can only take into consideration one of the two clusterings (repeated measures vs workplace) which one would you recommend q3. I have about 10000 observations and about 400 of them have missing values for the workplace. What would you recommend to do in this case? (also i should mention that the 205 workplaces are nested in 17 Organizations - For the moment i use general categories based on the organization: e.g. Organization1-Unclassified). Is there a meaningful way to actually impute these categories? q4. Would you recommend to use all 180 items for imputation or the items that i intend to use in each of my models? I use R for analysis and it would be greatly appreciated if you can recommend any packages for multiple imputation for clustered data. Thanks in advance


title: Dynamic estimation of required sample size in R

text: I have a prospective study with no data about estimated results, that could be used to get required sample size. Data looks like this:               caseID;groupID;value,result     1;1;12.3;0     2;1;15.6;1     3;2;11.3;0     4;2;13.4;1     ...      Is it possible to determine how much observations should be made to complete the study after adding new portion of data? I need both descriptive statistics ( % of cases have sign X, _result_ in example ) and comparative (2-4 groups compared by continuous variable, _value_ in example). Probability of type I error ($\alpha$) is 10% and type II ($\beta$) - 5% ( _for example_ ) Statistical significance and power - 0,9 and 0,95. The best solution is to have a grid like this for making desision to stop or continue the study:               a;b;n     5;5;30     10;5;28     5;10;29     ...      Where a is $\alpha$ , b is $\beta$ and n is the size of sample required for study with this $\alpha$ and $\beta$ (Type I and type 2 error probability). My data is imported to R. Is it possible to perform some calculations to get?:   1. Number observations required to complete study   2. Required number of observations for each group So I need to estimate sample size based on actual data to know it is time to finish the study. The question is: **How to do this in R, if it is possible?** Any suggestions are welcomed.


title: Quantile transformation in R

text: Is there a way to easily create factors based on quantiles of selected variables in a dataframe? Say, in a datatset D, I have variables V1 to V10, which are all numeric. I would like to create dummies for V7 to V10 based on their respective quantiles.


title: Stepwise regression vs. elastic net

text: I understand that Stepwise regression analysis has lots of limitations, including the assumption that the predictors are not highly correlated with each other. In fact, this limitation was the most important reason that I switched to Elastic Net, as I had 75 predictors in my model, some of which are highly correlated. Using Elastic Net, I could reduce my predictors to 21. I used these selected 21 variables in a multilinear regression model and calculated the coefficient of determination ($R^2=0.58$). However, when I used Stepwise analysis on the same data, only 11 variables got selected, while the R-square stayed the same! Does it mean that my results from Stepwise analysis can explain a higher proportion of my outcome? If so, how can I justify the limitations of Stepwise analysis over Elastic Net when I'm getting better results?


title: Testing marginal homogeniety with multivariate outcomes and matched pairs

text: Suppose I have 20 people matched in 10 pairs. Suppose two doctors are assigned at random to diagnose one person in each pair, so in each pair each person is diagnosed by a different doctor. Suppose that each doctor reports a vector of symptoms for each patient. That is, the response for patient $i$ diagnosed by doctor $j$ is a report $x_{ij}=(a,b,c,)^T$ where each element reports the number of instances of condition a, b, or c in patient $i$ like $x_{ij}=(2,5,0)$. Since the doctors are assigned to patients at random we would expect the average report $x_{.j}=(\bar{a},\bar{b},\bar{c})^T$ across doctors to be identical _if_ both doctors have the same diagnostic ability. How do I test the null that both doctors have the same diagnostic ability against the alternative that they differ in some respect?


title: Statistical test for two-factor heteregeneously distributed samples

text: **Data** I have a sample size of 50 in group A and 50 in group B where groups A and B are unmatched. Each sample in group A and group B has two frequencies associated with it, which I'll call $x$ and $y$. A histogram of $x_{A}$ shows that it is clearly Gaussian and a Shapiro test confirms this. However, the distribution of $x_{B}$ is - according to a Shapiro test - not Gaussian (p<0.01). For $y$, the distributions of group A and of group B are clearly not Gaussian but both have more of a Gamma-like or Weibull distribution (not sure which or whether it matters). **Objective** I want to be able to say whether any difference in $x_{A}$ and $x_{B}$ (as determined by some test statistic) is more significant (lower p-value) than any difference between $y_{A}$ and $y_{B}$. **Thoughts** Ideally, in order to make the p value from test with $x$ and the p value from test with $y$ comparable, I would think that the same statistical test would be needed. I know that when the spreads/distributions of $x_{A}$ and $x_{B}$ (and similarly for $y$ are the same, I can use the Mann Whitney test. However, the spreads (as described above are not the same). Therefore, could I do some sort of data transformation and do a t-test, or should I use a 2-sample Kolmogorov Smirnov test? In the latter case, how can I extend it to the _k_ -sample case?


title: Using standardized Y in Elastic Net

text: I have an Elastic Net model that is selecting a number of variables from X, for prediction of Y. The assumption for Elastic Net is that X is standardized (I'm using Z-Scores), and Y is centered around zero (I'm using Y-mean(Y)). So, I am wondering if my Elastic Net model will act differently if I use both standardized X and Y (i.e., z-scores for both X and Y)?


title: Strict Bonferroni - looking for p value

text: I run 18 exploratory regressions and would like to know what the strict Bonferroni p-value would be. Nine regressions are negative binomial regressions with 4 predictor variables, the other 9 are WLS regressions with 6 predictor variables. These associations are yet unexplored (my dependent variables are clinical questionnaire items instead of clinical sum scores, we are looking for heterogeneity). I am sure more information is needed, and would be glad if you could either help me directly (I will gladly provide all information that is necessary) or hint me at literature that enables a non-stats person to do the calculation himself. I know that Bonferroni is very strict and do not plan to actually apply it to my analyses - I plan to use .01 as type I error to detect significant effects (and would appreciate feedback). My dependent variables are substantially correlated (between .3 and .6), but I at least want to provide what the p-value would be in the paper and discuss why it would be too conservative.


title: Tail inequality on sum of product of normal variables

text: For independent random variables $ x_1,..,x_n$ and $y_1,...,y_n$ following normal distribution $N(0,1)$, I need a simple estimate formula for $P(| \sum_1^n x_iy_i | \leq nt ) \leq e^{(?)}$ for $t>1$. Thanks.


title: Intuition behind perceptron algorithm with offset

text: I was looking for an intuition for the perceptron algorithm with offset rule, why the update rule is as follows: _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)}\\\ \quad\theta^{(k+1)}_0 = \theta^{k}_0 + y^{(t)}\\\ \\}$ When the offset is zero, I think the update rule is completely intuitive. However, without it, it seems a little odd just adding 1 or -1 to the offset. The only reason I could come up with to explain it was the following but I don't really think its very intuitive explanation and was looking for a different explanation. My non-intuitive answer: When the perceptron makes a mistake then: $y^{(t)}(\theta^{T}x + \theta_0) \leq 0$ But we can re-write the top part as: $<\theta, \theta_0> \cdot <x^{(t)}, 1> = \theta^{T}x + \theta_0$ and now if we just appeal to the original perceptron rule and change the feature vector to have the one attached at the end and the normal now includes $\theta_0$, now the update would occur as following: $ \theta'^{(k+1)} = \theta'^{(k)} + y^{(t)}x'^{(t)}$ which is: $<\theta, \theta_0> \+ y^{(t)}<x^{(t)}, 1> = <\theta + y^{(t)}x^{(t)}, \theta_0+y^{(t)}>$ I think this might be correct, but even if it is, I didn't really think it was intuitive or "obvious" and was wondering if anyone had a different argument? Thanks! PS: Feel free to edit my algorithm to have indentation and spaces, I couldn't make it have indentation without losing the latex :(


title: Is it meaningful to use Kolmogorov-Smirnov test in small sub sample?

text: I have a data set with around 14,000 observations. This is a sample data set not the population. I have fitted a model (based on fuzzy logic) based on this data. I have 27 actual values(from 14,000 observations) ans 27 values by model. Is it meaningful to use Kolmogorov-Smirnov test to justify my model values and 27 actual values has come from a same distribution?


title: Expectation of Fisher information of log likelihood function

text: We evaluate Fisher information of a log likelihood function for a particular parameter in this process. How can we take expectation of the second derivative of log likelihood function with respect to a certain parameter? Which function is to be used while evaluating expectation?


title: Can p and q in ARMA(p,q) be estimated from its ACF and PACF?

text: From Wikipedia For AR(p), its p can be estimated from where its sample partial autocorrelation function (PACF) plot becomes zero. For MA(q), its q can be estimated from where its sample autocorrelation function (ACF) plot becomes zero. For ARMA(p,q), p > 0 and q > 0, its ACF and PACF are both decaying but not to zero. Can its p and q be estimated from its ACF and PACF? Thanks!


title: Reconstruct a "blocky" picture?

text: Consider a finite set $A$. Let the sample space be $A\times A$. We have an unknown probability distribution $f$ on this sample space. Now this probability distribution has a "blocky" property, which I am going to explain. Let a block be a function (not necessarily probability distribution) over the sample space $b_{(B,C,r)}$ where $B,C\subset A$ and $r$ be any real number (not necessarily in $[0,1]$). Then $b_{(B,C,r)}(x,y)=r$ if $(x,y)\in B\times C$ and $b_{(B,C,r)}(x,y)=0$ if not. Then $f$ can be written as a finite sum of those blocks (at most $|A|^{2}$ of them in fact). Let $X$ be the minimum number of blocks needed to sum to $f$. If prior to observing any trials, we know that $X$ follows a geometric distribution with unknown $p$ truncated to $[1,|A|^{2}]$ then we say that $f$ is "blocky". (Note that despite the analogy to reconstructing a picture, a block need not be a rectangle: there is no relationship among elements in $A$) Now we are given the outcome of a number of trials, where each trial consists of picking an element in the sample space with probability given by that unknown but "blocky" probability distribution $f$. The task is to make estimate of $f$ of each possible outcome. So what is a good algorithm for that? I'm sorry that I am new to statistics and thus I was not able to phrase this in a more proper way, and I do not even know if this problem has a solution or not, so it might be trivial or unsolved for all I know. Hope I can get some help on this. Thank you for any help you can give.


title: Scoring validation database based on estimates from test database in R

text: I divided my dataset into Test and Validation (50-50 split). I ran glm function (link=binomial) on Test dataset and got the parameter estimates. How do I score the Validation dataset based on these parameter estimates (beta) that I got from Test dataset. I know it has something to do with apply () but I am not sure. please advise.


title: About Standardizing Variables

text: Before data analysis, by $\frac{(X-mean)}{std}$, so the new variables have mean 0 and variance 1, and then compare different variables since then they are on the same measure scale, my question is, what is preserved about the data after standardizing then? Thanks!


title: Predicting with logit coefficients

text: I have a very beginner question about Logit coefficients. I am trying to predict how often the home team with win a basketball game. I did a regression analysis much like the one on this page and below are the variables and their coefficients. I haven't done a logistic regression and I'm sure how to utilize the coefficients to make predictions.               home_cumORTG     0.0766           away_cumORTG    -0.0891           home_cumDRTG    -0.1029           away_cumDRTG     0.0690           intercept        5.7061           Assume the home and the away team both have offensive and defensive ratings of 105 what is the probability that the home team will win? This is my best guess on how to do it but I hope it's not right because if it is my data is no good.               odds = 105 * .0766 + 105 * -.0891 + 105 * -.1029 + 105* .0690 + 5.7061      odds / (1+odds)           output 0.45477345837195365      Something is wrong because it is predicting that if two evenly matched teams play the home team will win less than half of the time when in reality the home team wins over 60 percent of the time. EDIT I think this is how it should be               1/(1+(1/ math.exp(5.7061 + 105 * .0766 + 105 * -.0891 + 105 * -.1029 + 105* .0690)))


title: Why does rejecting the null in goodness-of-fit tests not imply accepting the null?

text: From All of Statistics by Wasserman: > **Goodness-of-fit testing** has some serious limitations. If reject $H_0$ > then we conclude we should not use the model. But if we do not reject $H_0$ > we cannot conclude that the model is correct. We may have failed to reject > simply because the test did not have enough power. This is why it is better > to use **nonparametric methods** whenever possible rather than relying on > parametric assumptions. My questions are:   1. Are goodness-of-fit tests parametric or nonparametric?   2. why "if we do not reject $H_0$ we cannot conclude that the model is correct"?  If it is because "We may have failed to reject simply because the test did not have enough power", why "it is better to use nonparametric methods whenever possible rather than relying on parametric assumptions"? Doesn't the same reason and conclusion apply to nonparametric methods?   3. Is it correct that "Goodness-of-fit testing" here mean testing if the distribution of a sample is a specific distribution?  Does the same conclusion "if we do not reject $H_0$ we cannot conclude that the model is correct" also apply to testing if the distribution of a sample and the distribution of another sample are the same, such as z test for two normally distributed groups of samples, and Kolmogorov-Smirnov two-sample test? Thanks and regards!


title: Estimating a joint distribution from observed max and min samples

text: Suppose that you have jointly distributed $N$ (~100) random variables, $\\{X_1,\ldots,X_N\\}$, and this distribution is unknown to you. However you do know that their sum is zero by construction. Having $L$ (~3000) observations of $\max(X_1,\ldots,X_N)$ and $\min(X_1,\ldots,X_N)$ each, how can you make a statement about the joint distribution of $\\{X_1,\ldots,X_N\\}$? I gave it a shot by assuming that $\\{X_1,\ldots,X_N\\}$ have a joint normal distribution but when I looked at the QQ-plot of the observed samples and the simulated ones, I wasn't really convinced that the underlying distribution was normal. Even though I have no evidence that supports this hypothesis, I think I can assume that the marginal distributions of $X_1,\ldots,X_{N}$ are identical. From this assumption it follows that the mean of $X_i$ is zero for all i running from 1 to N. Furthermore, for the covariance matrix we can write   $\sigma_{ij} = \sigma^2$ whenever $i = j$   $\sigma_{ij} = -\frac{\sigma^2}{N-1}$ whenever $i \neq j$ Here $\sigma_{ij}$ is the $ij$th element of the covariance matrix and $\sigma^2$ is the variance of $X_i$ for all $i$ running from 1 to N.


title: Testing prediction time series against real data

text: Say I have a series of forecasts and observations like this:                     EntityF   EntityO     2004  120       125     2006  166       173     2008  150       167        2010  152       -      And assume that the (i) entity is the same and (ii) the forecasting methodology is constant. I'd like to   1. Produce a meaningful metric of the forecasting error.   2. Be able to predict the current forecast (2010) error based on 1.


title: One-way ANOVA with unequal sample sizes

text: I have a function for performing one-way ANOVA in R:               cond<-gl(4,5,20,label=c("a","b","c","d"))          aof<-function(x){     m<-data.frame(cond,x);     anova(aov(x~cond,m))     }          anova.results<-apply(x,1,aof)      I use it to perform ANOVA (testing four conditions, each with 5 samples) on several hundred rows of `x` (genes). However, I need to drop some samples such that two of my conditions have 5 samples and two have 4. I'm having trouble reworking my function to allow for ANOVA with unequal sample sizes. Any help is greatly appreciated! Edit for clarification: I am using R, but no special package, just the standard `stats` package. From what I understand, the `gl` function specifies the factor level pattern. `gl(n,k,n*k,labels)` such that `n` is the number of levels (four conditions, in my case), `k` is the number of replications (5 for my original analysis), and `n*k` is the total number of observations. `x` is a matrix with 20 columns, such that the five samples for condition `a` are first, then condition `b`, etc. The rows have values for each sample, and the `anova` function is done to each row individually. I have tried changing `k` to `c(5,4,4,5)` for unequal sample sizes, but R says it is an `invalid 'times' value`. So, I'm trying to figure out another way to rework this function, or write another ANOVA function.


title: Binom.test error in R: " formal argument "p" matched by multiple actual arguments"

text: I'm running a binom.test on the data set UCBAdmissions (comes with R) and am stuck on an error message. **About the data:**               > str(UCBAdmissions)     table [1:2, 1:2, 1:6] 512 313 89 19 353 207 17 8 120 205 ...     - attr(*, "dimnames")=List of 3     ..$ Admit : chr [1:2] "Admitted" "Rejected"         ..$ Gender: chr [1:2] "Male" "Female"     ..$ Dept  : chr [1:6] "A" "B" "C" "D" ...      *Number of females: 1835 *Number of femaled Admitted: 557 *Number of males: 2691 *Number of males Admitted: 1198 **I'm trying to answer this question:** "Do you think there was gender discrimination?" **It's a Bernoulli distribution so I am using the binom.test, however I get an error message that I can't make sense of:**               no.w <- sum(UCBAdmissions[1,2,])     no.acc.w <- sum(UCBAdmissions[1,2,])     binom.test(no.acc.w,no.w, p = .304, c="two-sided", p=.95)     Error in binom.test(no.acc.w, no.w, p = 0.304, c = "two-sided", p = 0.95) :      formal argument "p" matched by multiple actual arguments      In advance, thanks for poking holes in this and pointing me in the right direction. Jen


title: How to correct for non-linearity of response in linear regression

text: I want to train a linear regression model to predict a non-linear variable. This how the two independent variables correlated against the response (points are jittered): ![enter image description here](http://i.stack.imgur.com/U1Ksu.jpg) ![enter image description here](http://i.stack.imgur.com/Hf95A.jpg) And the residuals against the fitted values: ![enter image description here](http://i.stack.imgur.com/itQfs.jpg) Most of the values for the response are zero. The effect is a very strong heteroscedasticity                       studentized Breusch-Pagan test          data:  model     BP = 55483.84, df = 2, p-value < 2.2e-16      event though the the predictors are strongly correlated with the response               Call:     lm(formula = response ~ predictor1 + predictor2, data = train_predictors)          Residuals:         Min      1Q  Median      3Q     Max      -7.6996 -0.0268 -0.0238 -0.0182  4.8785           Coefficients:                   Estimate Std. Error t value Pr(>|t|)         (Intercept)  2.748e-02  2.825e-04   97.28   <2e-16 ***     predictor1   8.491e-05  6.574e-07  129.16   <2e-16 ***     predictor2  -3.934e-10  8.298e-12  -47.41   <2e-16 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1          Residual standard error: 0.1561 on 498498 degrees of freedom     Multiple R-squared:  0.0365,    Adjusted R-squared:  0.0365      F-statistic:  9442 on 2 and 498498 DF,  p-value: < 2.2e-16      Should I consider more adopting non-linear models or could I first try correcting the non-linearity of the response?


title: Understanding random samples

text: I was looking at the following statement and felt I needed some help in sorting it out: A random sample of size n from a population $f(x)$ is a collection of n independent random variables $X_1,...,X_n$, each having the distribution $f(x)$. let us apply the above statement to all freshmen at a certain university where we are interested in their heights and weights, So I guess:   1. $f(x)$ would be the population of freshmen at a certain university.   2. What exactly does $x$ signify in $f(x)$   3. so $X_1,...,X_n$ would be the weight measured for a sample of size $n$ ?    4. so for heights would be need $Y_1,...,Y_n$?   5. strictly speaking can we use $X$ for both height and weight (I guess no -- it would be meaningless as we're mixing types)?   6. The last part of the statement "each having the distribution $f(x)$" has me a bit on tender ground -- would be interested in getting clarifications about that. Note: The statement above is from: "Statistical Concepts and Methods (Wiley Series in Probability and Statistics), Gouri K. Bhattacharyya, Richard A. Johnson", page 208.


title: Regression discontinuity versus matching with spatial discontinuity

text: I am interested in using a spatial research design. Imagine a line, like a time zone line. For example, in the United States, the line that makes between Eastern Standard Time and Central Standard Time runs North to South through the U.S. (and other places), more or less. Suppose the United States implemented a policy where, with respect to county boundaries, all counties west of the 75th longitude were to get a tax deduction. No counties to the right of it would. I want to evaluate the impacts of the this tax deduction on crime. For the sake of argument, can we please assume that the discontinuity is valid. I am trying to understand this idea for spatial discontinuities, in general? I can imagine proceeding one of two ways. One way would a regression discontinuity. Exploiting the distance of say a county from this line. The other way would be a matching of counties that straddle this line on opposite sides. Note that a county might appear more than once if it borders more than one county on the opposite side of the time zone line. The regression discontinuity includes a larger sample because it may include counties that are not exclusively along the border. The matching hones in on counties that touch each other, whereas the regression discontinuity, at best, can partition the time zone into smaller segmented lines where segments likely include multiple pairs of counties along the time zone. What are the trade-offs of the research design? Which is preferred? Does one offer something that the other doesn't? Thank you very much for your help!


title: How to do factor analysis for categorical data?

text: Now I have two data set which is generated by two listing flow of a e-commerce site. And these data set has different conversion rate, my task is to analysis the two data set and found out which factors lead the different conversion rate. BTW, most of the data is categorical data. Could you help point me some tutorial or method that I can refer to do the analysis ? Currently what I did is just according my intuition, I'd like to find out to some systematic method for factory analysis of categorical data.


title: Auto-correlation of Random numbers

text: Is saying random numbers are independent equivalent to saying random numbers are auto-correlated in simulation. I am using Auto-correlation method in simulation.


title: KL divergence between 2 distributions with unequal cardinalities?

text: Say $X$ is a discrete random variable with cardinality $|X|$ and $Y$ is a discrete random variable with cardinality $|Y|$. Does it make sense to talk about the KL divergences $D_{KL}(X||Y)$ or $D_{KL}(Y||X)$ of these 2 probability distributions if $|X| \neq |Y|$ ? If so, how does one compute it ?


title: Once again: paired versus unpaired t-tests

text: This is much the same question that I asked a few weeks ago, but I hope to explain myself more clearly this time. I start with 40 mice. I use my own scheme, based on my own ideas for how to make the two mice in a pair as similar as possible in terms of the experiment I am about to do. Pairing the 40 mice gives me 20 pairs; from now on the pair structure is unchanged. Using a random number generator, I choose one mouse from each pair, selected to be given the drug A; the other mouse in each pair is given drug B. From now on, the two mice in a pair are treated as similarly as possible in all respects, except that they are given the two different drugs A and B. For example they are always fed at the same time. We also take the precaution of doing the experiment in a blinded fashion, that is, arranging the experiment so that only one person X knows which mice are getting A and which B, and so that X has no other role in the experiment. The experiment has a numerical outcome, namely the cholesterol concentration in the blood at the end of the experiment, so we get 20 numbers. Let's assume that two distributions involved are approximately normal, with the same variance. The null hypothesis is that the means are the same. The mean for mice given drug A mice turns out to be lower than the mean for mice given drug B. Moreover, a two-sided paired t-test gives me a p-value of 0.005, which is eminently publishable, and an two-sided unpaired t-test gives a p-value of 0.06, which would make the work unpublishable according to criteria that are standard in biology publications. An eminent expert in the application field says that my pairing scheme is "not biological", and that the "correct" p-value should therefore be 0.06. My own reaction is:   1. that it is possible that the result is a fluke, so maybe the experiment should be repeated, particularly if the conclusions seem unlikely to experts.   2. Barring statistical flukes, what has been shown is that drug A is more effective at lowering the blood cholesterol in a mouse than drug B, and this result is statistically significant.    3. Barring flukes, the experiment shows that my very own pairing method has a sound biological basis, though an understanding of why the method is sound may still be unavailable.   4. Unless something is wrong with the experiment, other than the pairing scheme, the eminent expert is wrong.  What is the response of the Cross-Validated community? I'm not asking for practical advice about what to do with the experimental results---the paper is already published. I'm wanting to make quite sure that what I say is, in principle, correct, as it will affect my advice for future experiments.


title: Frequency of a wavelet filter

text: 1. Assume that I have decomposed a data set using Symlet Wavelet with six levels. How can I estimate the approximate frequency interval of each level? That would be great if you consider your answer in the Mathematica environment.   2. If the above dataset describes turbulent velocity field, how I can decide up to which level is the bulk velocity? I need to separate the turbulent fluctuations from the mean. The sampling frequency is 1 Hz and the duration of observation is 25 hours.


title: Testing if data comes from a normal distribution with mean 0 and unknown variance in Matlab

text: ## Question I have a vector of data, and I want to test if it came from a normal distribution with mean zero and unknown variance. Do you know if there is matlab function or simple script for this? If you don't know anything matlab specific, then a name and reference for the specific test is fine and I will just implement it myself. Also, if the specific test can return the confidence level instead of just answering yes-no at a given confidence level then that would be a benefit, but is not essential. * * * ## What I already know If I want to test if my data is from a normal distribution with mean 0 and variance 1 then I can use the Kolmogorov-Smirnov test. If I want if my data is from a normal distribution with unknown mean AND variance then I can use the Lilliefors test or the Jarque-Bera test. However, I want a fixed mean (= 0) and unknown variance. ## Naive approach The naive approach is to take my data $D$, calculate the variance around zero $\sigma^2_0$ and then renormalize my data by this to get a dataset $D'$. Then I can perform the Kolmogorov-Smirnov test on this. However, it is not clear how one would justify this, especially since the KS tests specifically warns against testing against distributions with parameters estimated from the same data (renormalizing $D$ to $D'$ will be the same as testing against a normal distibution with mean zero and variance $\sigma^2_0$). Is this naive approach justified?


title: How to test for differences from baseline where dependent variable is non-normally distributed using SPSS?

text: I have some physiological measurements which consist of 4 events measured in a "during task" condition and preceded and followed by a rest period. I want to test that the data points immediately preceding the 4 events do not differ from the rest baseline, calculated as the average of the values in the start and end rest periods. My problem is that I have measurements from 13 subjects, 3 session each, each session consisting of 9 of these measurement blocks. To make all the data comparable, I subtracted the baseline from the values preceding each event, so that the test would now be against zero for all of them. Given that the distribution of the values is not normal, how should I conduct a test in SPSS (PASW v. 18) to prove that overall these data points are not significantly different from zero, resulting in one statistic to report?


title: What are good ways of visualizing many effects ordered in groups and subgroups?

text: In epidemiology, we often deal with lots of factors associated with disease. This multitude of factors makes plots of effect sizes (often hazard ratios) of individual factors confusing. Accepting some oversimplification, factors can for the most part be reasonably assigned to some group like nutritional, lifestyle, social, psychological, demographic, and maybe even subgroups like e.g. macro-/micronutrients for nutritional factors. I am searching for a plot that gives an overview of the effects within these groups while facilitating a visual comparison of size and distribution of effects within and between groups. I have tried   * small multiples of histograms/smoothed density estimates and parallel boxplots of HRs in each group   * forest-plot-like plots showing effects sizes and their CIs as points and horizontal segments and groups separated by horizontal lines These ad-hoc approaches are not optimal. The first discards most of the information on the level of the individual factor, the second is too detailed on that level. Neither provides real subgrouping. Neither provides an informative visual comparison of distribution AND size of effects within groups. Is there something better?


title: Cobb-Douglas production: interpret regression

text: If we take the logarithm of the Cobb-Douglas production function, we get: ln(Y)=A+$\beta_1$ln(L)+$\beta_2$ln(K)+$\epsilon$ln(e) I understand that in the production function, the coefficients $\beta_1$ and $\beta_2$ are output elasticities of labor and capital. What would be the interpretation of these coefficients in the regression - is it just the change in output when the labor or capital share changes Should the beta coefficients add up to 1, and if yes, are there any issues of them changing simultaneously? Thank you!


title: How to get both MSE and R2 from a sklearn GridSearchCV?

text: I can use a GridSearchCV on a pipeline and specify scoring to either be `'MSE'` or `'R2'`. I can then access `gridsearchcv._best_score` to recover the one I specified. How do I also get the other score for the solution found by GridSearchCV? If I run GridSearchCV again with the other scoring parameter, it might not find the same solution, and so the score it reports might not correspond to the same model as the one for which we have the first value. Maybe I can extract the parameters and supply them to a new pipeline, and then run `cross_val_score` with the new pipeline? Is there a better way? Thanks.


title: What are good references for Deming's regression?

text: Do you know some good references (papers or books) for the theoretical inference in Deming's regression model ? **EDIT:** I was a little disconcerted about a point in Ripley and Thompson's paper Regression techniques for the detection of analytical bias. I size the opportunity to share my clarifications (thanks to @Procrastinator for his help) The postulated model is $x_i = u_i + \delta_i$ and $y_i = v_i + \epsilon_i$ where $u_i$ and $v_i=\alpha+\beta u_i$ are fixed numbers and the error terms $\delta_i$ and $\epsilon_i$ are independent. This is the commonly called _functional relationship model_. See Cheng & Ness's paper "Structural and functional models revisited" for a thorough review of the functional relationship model. But this is also nothing but Deming's well-known regression model. After substituting $x_i-\delta_i$ for $u_i$ in the expression of $y_i$ one gets $y_i=\alpha+\beta x_i + (\epsilon_i - \beta \delta_i)$ and Ripley and Thompson consider $y_i=\alpha+\beta x_i + (\epsilon_i - \beta \delta_i)$ as a weighted least squares regression model with covariate $x_i$ and error term $\epsilon_i - \beta \delta_i$. The sentence on page 349 _"... showing dependence of $y_i$ on $x_i$"_ is quite disconcerting: $x_i$ and $y_i$ are independent in the original postulated model. In fact the authors do not consider that $y_i=\alpha+\beta x_i + (\epsilon_i - \beta \delta_i)$ provides a _valid_ weighted least squares model but they do not explain _why_ this model is not valid. The cause, which is not clearly mentioned in the paper, is that $x_i$ is not independent of the error term $\epsilon_i - \beta \delta_i$ (I guess there are more details in Ripley's unplublished paper mentioned in the references). In fact one of the purposes of this paper is to discuss about the validity of "approximately valid" models. I'm still looking for a reference with a thorough presentation of the functional model and related models (structural and ultrastructural). Cheng & Ness's paper "Structural and functional models revisited" is a nice rigorous review about the theoretical treatment of the model but it does not address the more practical questions.