title: Describing the non-summative combination of proportions

text: I'm sure this an obvious question for someone. When one combines proportions, but the overall proportion cannot exceed 100 % in reality - what is this called? For example, if I'm trying to reduce the likelihood of someone catching a disease and a medication reduces the probability by 80% and following a meal plan reduces it by 40% (independent studies no interaction measured) then the overall reduction is worked out as follows: 80% + (40% of 20%) = 88%. I may have reduced the likelihood of someone catching the disease by 88%. How do I explain this mathematically? I know there's a probabilistic basis to this (is it the independent events side of things in basic probability, is it something to do with asymptotes, is it something to do with the math behind ANOVAs (i.e. interactions and independent)?). Is it all of these things combined?! I'm happy to head off on my own for further research I just need to know I heading to the right literature as i'm not a natural statistician! Any help is appreciated.


title: Determining sample size necessary for bootstrap method / Proposed Method

text: I know this is a rather hot topic where no one really can give a simple answer for. Nevertheless I am wondering if the following approach couldn’t be useful. The bootstrap method is only useful if your sample follows more or less (read exactly) the same distribution as the original population. In order to be certain this is the case you need to make your sample size large enough. But what is large enough? If my premise is correct you have the same problem when using the central limit theorem to determine the population mean. Only when your sample size is large enough you can be certain that the population of your sample means is normally distributed (around the population mean). In other words, your samples need to represent your population (distribution) well enough. But again, what is large enough? In my case (administrative processes: time needed to finish a demand vs amount of demands) I have a population with a multi-modal distribution (all the demands that are finished in 2011) of which I am 99% certain that it is even less normally distributed than the population (all the demands that are finished between present day and a day in the past, ideally this timespan is as small as possible) I want to research. My 2011 population exists out of enough units to make $x$ samples of a sample size $n$. I choose a value of $x$, suppose $10$ ($x=10$). Now I use trial and error to determine a good sample size. I take an $n=50$, and see if my sample mean population is normally distributed by using Kolmogorov-Smirnov. If so I repeat the same steps but with a sample size of $40$, if not repeat with a sample size of $60$ (etc.). After a while I conclude that $n=45$ is the absolute minimum sample size to get a more or less good representation of my 2011 population. Since I know my population of interest (all the demands that are finished between present day and a day in the past) has less variance I can safely use a sample size of $n=45$ to bootstrap. (Indirectly, the $n=45$ determines the size of my timespan: time needed to finish $45$ demands.) This is, in short, my idea. But since I am not a statistician but an engineer whose statistics lessons took place in the days of yonder I cannot exclude the possibility I just generated a lot of rubbish :-). What do you guys think? If my premise makes sense, do I need to chose an $x$ larger than $10$, or smaller? Depending on your answers (do I need to feel embarrassed or not? :-) I'll be posting some more discussion ideas. response on first answer Thanks for replying, Your answer was very usefull to me especially the book links.   But I am affraid that in my attempt to give information I completely clouded my question. I know that the bootstrap samples take over the distribution of the population sample. I follow you completely but... Your original population sample needs to be large enough to be moderately certain that the distribution of your population sample corresponds (equals) with the 'real' distribution of the population. This is merely an idea on how to determine how large your original sample size needs to be in order to be reasonably certain that the sample distribution corresponds with the population distribution. Suppose you have a bimodal population distribution and one top is a lot larger than the other one. If your sample size is 5 the chance is large that all 5 units have a value very close to the large top (chance to ad randomly draw a unit there is the largest). In this case your sample distribution wil look monomodal. With a sample size of a hundred the chance that your sample distribution is also bimodal is a lot larger!! The trouble with bootstrapping is that you only have one sample (and you build further on that sample). If the sample distribution really does not correspond with the population distribution you are in trouble. This is just an idea to make the chance of having 'a bad sample distribution' as low as possible without having to make your sample size infinitly large.


title: Proving significantly better performance in binary forecasting

text: I have two different algorithms that make forecasts for binary events. The observed result can either be 1 or 0 (like "rain" or "no rain"). The algorithms usually give a forecast in the 0.4-0.6 range. For measuring error I take the difference between forecast and actual result to the power of 2. So if algorithm 1 said there was a 45% chance of raining and it did rain the resulting error from that observation would be (1-0.45)^2. I do this for all observations and sum up the errors over all observations (giving me something like the Brier Score, with the only difference being that I'm not dividing by #of observations). Suppose algorithm 1 has a lower sum of error than algorithm 2. How can I test whether algorithm 1 is _significantly better_ than algorithm 2?


title: What is the standard deviation of a gaussian filter?

text: I was referring to this thesis, High dimensional gaussian filtering for computational photography (or pdf), for understanding gaussian filters. However, I didn't understand one part. The gaussian filter is given by $$ \hat{v_i} = \sum_{j}\exp(-|p_i-p_j|^2 / 2)\,v_j $$ If $p_i$ and $p_j$ are position vectors then it is simply the Gaussian blur. Now I didn't understand if we divide $p_i$ and $p_j$ by a certain constant, the standard deviation of the filter is changed. I didn't get how it occurs. Also I didn't get what you mean by the standard deviation of a filter?


title: How to collapse ordinal variable values into significant different values?

text: I have an ordinal variable related to an outcome that is comprised of many levels and I´d like to collapse the number of ordinal values as much as possible.               > require(ipred)     > require(party)     > data(GBSG2)     > head(GBSG2)       horTh age menostat tsize tgrade pnodes progrec estrec time cens     1    no  70     Post    21     II      3      48     66 1814    1     2   yes  56     Post    12     II      7      61     77 2018    1     3   yes  58     Post    35     II      9      52    271  712    1     4   yes  59     Post    17     II      4      60     29 1807    1     5    no  73     Post    35     II      1      26     65  772    1     6    no  32      Pre    57    III     24       0     13  448    1     > table(GBSG2$tgrade)            I  II III       81 444 161      > ctree(Surv(time,cens)~tgrade,data=GBSG2) -> mn     > plot(mn)      ![enter image description here](http://i.stack.imgur.com/WYIUd.png) Would it be correct to claim that `tgrade` here could be collapsed into two instead of three values? edit: Running the usual parametric analysis I get:               >     anova(i1,i2)     Analysis of Deviance Table      Cox model: response is  Surv(time, cens)      Model 1: ~ tgrade      Model 2: ~ tgrade == "I"        loglik  Chisq Df P(>|Chi|)       1 -1776.0                           2 -1778.1 4.3049  1     0.038 *     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1      >     anova(i1,i3)     Analysis of Deviance Table      Cox model: response is  Surv(time, cens)      Model 1: ~ tgrade      Model 2: ~ tgrade != "III"       loglik  Chisq Df P(>|Chi|)         1  -1776                             2  -1784 16.033  1 6.225e-05 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1      > anova(i2,i3)     Analysis of Deviance Table      Cox model: response is  Surv(time, cens)      Model 1: ~ tgrade == "I"      Model 2: ~ tgrade != "III"        loglik  Chisq Df P(>|Chi|)         1 -1778.1                             2 -1784.0 11.728  0 < 2.2e-16 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1      > extractAIC(i1)     [1]    2.000 3555.975     > extractAIC(i2)     [1]    1.00 3558.28     >   extractAIC(i3)     [1]    1.000 3570.008        Hence the i1 model provides a better fit than i2 and i3, and i2 fits significantly better than i3. So now all three categories are warranted with respect to survival, which is at odds with the ctree approach. Can anyone explain this? Is this due to the conditional nature of ctree instead of the semiparametric nature of cox regression?


title: What is the difference between Informative (IVM) and Relevance (RVM) vector machines

text: I'm trying to understand if there is any specific difference between Informative IVMs and Relevance RVMs other than the terminology. I've not seen anything explicit. When I'm reading about vector machines, it is easy to see the difference of IVM/RVMs from Support (SVM) vector machines [colloquially, for classification, the SVM finds those points (vectors) that define the DMZ (de-militatized zone ;-) between the categories, while the RVM finds those that are the 'middle' of the crowd, and an associated crowd 'size' (e.g. in gaussian globs)], but I don't see any special difference between I/R vector machines beyond a choice of terminology by their proponents. Is there a difference?


title: Whether to use correlation matrix between latent variables or raw data as input to SEM?

text: I have a fairly large model - a total of 8 latent variables and each one has 10 items (observed variables). Drawing this model is a mess (there is a mediation effect as well). Is it better to just compute correlations among my latent variables and import the correlation table instead of using all the raw data?


title: Pointwise mutual information for tag generation in Twitter

text: I came across the following link. http://stackoverflow.com/questions/2764116/tag-generation-from-a-small-text- content-such-as-tweets The question is about how do we assign tag to tweets automatically, the answer suggest using pointwise mutual information, however, my question is if how can we achieve it? If we have 1000 unique tweets and 100 unique tags, the probability of hashtag $x$ appear in tweet $y$ $p(x,y) = 1$, $p(x) = 1/1000$ and $p(y) = 1/1000$ will always stay, right? Can someone please enlighten me on this please? Thanks. Regards, Andy.


title: AR process invertibility and stationarity conditions

text: I am studying the stationarity and invertibility conditions of AR processes. Concerning the stationarity condition of the simplest form of an AR(1) process (without constant, drift). Let us say that we convert such an AR process to an infinite MA process. If we take the expected value of the new form of the process (the infinite MA) then it is going to be zero since the error term is zero mean white noise so its mean value is 0. So an infinite sum of zeros will be zero (maybe its here that i am wrong). Also its variance is finite (sigma squared) so if we take the variance of the MA process it will be sigma squared divided by one minus phi squared(as an infinite sum of terms of a geometric progress). So since we can prove that the process has constant mean and finite variance without any condition why Diebold says that in order for the AR(1) process to be covariance stationary the AR coefficient should be less than one in absolute value? To be more specific, Diebold says that the MA representation of for y is convergent if and only if phis is less than 1 in absolute value. Thus this is the condition for covariance stationarity in the simple AR(1) case (without constant, drift). The only though that i have is that the sum of infinite zeros cannot be defined so this is why we need the MA form of the AR process to be convergent. And if it is convergent (if phi less than one in absolute value) then the summation of zeros will be finite and hence the expected value of the series will be constant (0). Is this the answer or it is something else?


title: Generate sets of values with high correlation coefficient

text: Apology if this is too simple. I couldn't get the more advanced r-help group to respond. I am planning to characterize workloads by measuring the correlation coefficient of two sets of real values but before that I wish to generate two sets of sample values that have a high coefficient and a low coefficient. I want to plot both in the same graph so that I can see the highly correlated values' together(peaks and troughs). I use R and know about rseek. If there is any particular R book that could help my capcaity planning efforts I will buy it. Generate a random variable with a defined correlation to an existing variable is a tad too advanced for me at this time. Note : The two sets of values that I am about to plot are related because I am plotting CPU usage and a througput number. So if the no: of bytes increases the CPU usage may increase. Both are postitive values. So if the correlation is high both will either increase together or decrease together. Thanks.


title: How to simulate functional data?

text: I'm trying to test various functional data analysis approaches. Ideally, i'd like to test the panel of approaches i have on simulated functional data. I've tried to generate simulated FD using an approach based on a summing Gaussian noises (code below), but the resulting curves look much too rugged compared to the real thing. I was wondering whether somebody had a pointer to functions/ideas to generate more realistic looking simulated functional data. In particular, these should be smooth. I'm completely new to this field so any advice is welcomed.               library("MASS")     library("caTools")     VCM<-function(cont,theta=0.99){         Sigma<-matrix(rep(0,length(cont)^2),nrow=length(cont))         for(i in 1:nrow(Sigma)){             for (j in 1:ncol(Sigma)) Sigma[i,j]<-theta^(abs(cont[i]-cont[j]))         }         return(Sigma)     }               t1<-1:120     CVC<-runmean(cumsum(rnorm(length(t1))),k=10)     VMC<-VCM(cont=t1,theta=0.99)     sig<-runif(ncol(VMC))     VMC<-diag(sig)%*%VMC%*%diag(sig)     DTA<-mvrnorm(100,rep(0,ncol(VMC)),VMC)            DTA<-sweep(DTA,2,CVC)     DTA<-apply(DTA,2,runmean,k=5)     matplot(t(DTA),type="l",col=1,lty=1)


title: Visualizing longitudinal data with binary outcome

text: For longitudinal data with a numeric outcome, I can use spaghetti plots to visualize the data. For example something like this (taken from the UCLA Stats site):               tolerance<-read.table("http://www.ats.ucla.edu/stat/r/faq/tolpp.csv",sep=",", header=T)     head(tolerance, n=10)     interaction.plot(tolerance$time, tolerance$id, tolerance$tolerance,                       xlab="time", ylab="Tolerance",  legend=F)      ![enter image description here](http://i.stack.imgur.com/a2E5i.png) But what if my outcome is binary 0 or 1? For example, in the "ohio" data in R the binary "resp" Variable indicates the presence of a respiratory disease:               library(geepack)     ohio2 <- ohio[2049:2148,]     head(ohio2, n=12)           resp  id age smoke      2049    1 512  -2     1      2050    0 512  -1     1      2051    0 512   0     1      2052    0 512   1     1      2053    1 513  -2     1      2054    0 513  -1     1      2055    0 513   0     1      2056    1 513   1     1      2057    1 514  -2     1      2058    0 514  -1     1      2059    0 514   0     1      2060    1 514   1     1          interaction.plot(ohio2$age+9, ohio2$id, ohio2$resp,                       xlab="age", ylab="Wheeze status", legend=F)      ![enter image description here](http://i.stack.imgur.com/c2M0g.png) The spaghetti plot gives a nice figure, but is not very informative and does not tell me much. What would be a suitable way to visualize this kind of data? Maybe something that includes a probability-value on the y-axis?


title: A question about parameters of Gamma distribution in Bayesian econometrics

text: The Wikipedia article on the Gamma distribution, lists two different parameterisation methods, one of them frequently used in Bayesian econometrics with $\alpha>0$ and $\beta>0$, $\alpha$ is shape parameter, $\beta$ is rate parameter. $$X\sim \mathrm{Gamma}(\alpha,\beta).$$ In a Bayesian econometrics textbook written by Gary Koop, the precision paramether $\frac{1}{\sigma^2}=h$ follows a Gamma distribution, which is a **prior distribution** $$h\sim \mathrm{Gamma}(\underline{s}^{-2},\underline{\nu}),$$ where $\underline{s}^{-2}$ is mean and $\underline{\nu}$ is degrees of freedom according to his Appendix. Also $s^2$ is standard error with definition $$s^2=\frac{\sum(y_i-\hat{\beta}x_i)}{\nu}.$$ Thus for me, these two definition of the Gamma distribution are completely different, since the mean and variances will be different. If we follow the wikipedia definition, the mean will be $\alpha/\beta$, not $\underline{s}^{-2}$. I am highly confused here, would anyone help me to streighten the thoughts here?


title: Calculating probability of project tasks, and most likely ending date?

text: Say I have a project (you know the typical one you can represent using a Gantt Chart ) If I where to assign a particula probability to each of the tasks in the Gantt... say for example:   * Task1 (2 days, i am 75% sure of that)    * Task2 (4 days, i am 50% sure of that) for it to start, Task 1 most be finished    * Task3 (8 days, i am 25% sure of that) for it to start, Task 2 most be finished How can I calculate how likely it is that I will finish this 14 day project? ( _In case you are wondering why I am asking this, I am a software developer, and I am constantly creating project plans and I would like to learn more on this subject so that I can project estimation better_ )


title: What is acceptable error in prediction models?

text: I have developed a model for construction waste reduction, but I'm not sure how to validate it to serve the purpose. I have two sets of data: measured waste % from sites and predicted values obtained from model constructed. How do I tell that the model constructed is good for prediction purposes? The data is linearly related.


title: Detect circular patterns in point cloud data

text: For some volume reconstruction algorithm I'm working on, I need to detect an arbitrary number of circular patterns in 3d point data (coming from a LIDAR device). The patterns can be arbitrarily oriented in space, and be assumed to lie (although not perfectly) in thin 2d planes. Here is an example with two circles in the same plane (although remember this is a 3d space): ![enter image description here](http://i.stack.imgur.com/FZxAU.png) I tried many approaches.. the simplest (but the one working best so far) is clustering based on disjoint sets of the nearest neighbor graph. This works reasonably well when the patterns are far apart, but less so with circles like the ones in the example, really close to each other. I tried K-means, but it doesn't do well: I suspect the circular point arrangement might not be well suited for it. Plus I have the additional problem of not knowing in advance the value of K. I tried more complicated approaches, based on the detection of cycles in the nearest neighbor graph, but what I got was either too fragile or computationally expensive. I also read about a lot of related topics (Hough transform, etc) but nothing seems to apply perfectly in this specific context. Any idea or inspiration would be appreciated.


title: Regarding importance of variables after applying PCA

text: I have a data with 15 variables. I did PCA on it to come up with 3 factors all of which had 5 variables each that explained the corresponding variance in those factors. After completing the entire process, I ran Logistic regression with these factors being my independent variables. I have individual weights of these factors (for prediction) ready as well. Now what I want to know is the individual importance of each variable in doing the subsequent predictions? As in what role (weight) each variable plays in predicting the dependent variable? Thanks in advance.


title: How to calculate relative error of a random variable?

text: This maybe a naive question, but I found two ways to calculate the relative error, one is common and another is from the papers I read. Assume that the exact number is $n$, and the estimated number is $\tilde{n}$, the two ways to calculate the relative error are: 1) relative error = $\frac{n-\tilde{n}}{E[n]}$, which is the standard definition 2) relative error = $\frac{\sqrt{Var(\tilde{n})}}{n}$, which is the coefficient of variance (I read from a paper). Now I am confused which one I should use, is there any corresponding scenario that these two formulas should be used? PS: Please note that $n$ is a random variable.


title: Order statistics: Probability that a normal random variable is k-th out of n when ordered

text: My question is almost the same as that of Order statistics: probability random variable is k-th out of n when ordered with the exception that the underlying distribution from which $X_1$ is drawn, having CDF $F$, is $N(\mu_1,\sigma_1^2)$. And the distribution from which $X_2,X_3,\ldots,X_n$ are drawn, having CDF $G$, is $N(\mu_2,\sigma_2^2)$. If that is the case then how to calculate the probability that $X_1$ is the $k^{th}$ out of the $n$ random variables when they are sorted, in order. Even if the probability might not be expressible by any elementary formula for large $n$ could you please give hints as to how to numerically approximate. Also, in that case, could you please solve for the case of some small $k$ and $n$.


title: How to assess whether a set of cells all have a specified value in Excel

text: I have collected a bunch of statistical results in the form of "YES" and "NO" strings. Now I would like to have a summary cell displaying "YES" if all these cells equal "YES" (or are empty).


title: In relation to clinical trials, what is clinical reasoning in contrast to statistical reasoning?

text: I read _Clinical Trials, a Methodologic Perspective_ (S. Piantadosi) as I was suggested by one of you. According to the author: > A trialist must understand two different modes of thinking that support the > science-clinical and statistical. They both underlie the re-emergence of > therapeutics as a modern science. Each method of reasoning arose > independently and must be combined skillfully if they are to serve > therapeutic questions effectively. I cannot figure out what the author means by "clinical reasoning". Can you help me to understand that notion? Thank you in advance


title: LIBSVM overfitting

text: I trained two svms (LIBSVM) with 15451 samples after I did a 10-fold cross- validation and found the best parameter values for gamma and C (RBF kernel). In one svm I used just 1 feature and in the second an additional one (to see whether this additional is improving prediction). After CV I have am accuracy of 75 % (SVM with one feature) and 77 % (SVM with that additional one). After testing on another 15451 instances I have an accuracy of 70 % and 72 % respectively. I know that this is called overfitting but is it significant here, since it is only a difference of 5 %. What could I do to avoid overfitting? Is it even good to use just one or two features and a relatively big training set? Hope you can help me out.


title: Does a psychological scale need to be revalidated after reducing the number of response options from 9 to 4?

text: I have used an established 33-item empathy questionnaire. I have modified the scale from a 9 point to a 4 point scale in order to force an answer and reduce ambiguity. **Do I need to retest the scale properties after this reduction in the number of response options?**


title: Difference between estimation and learning

text: What is the difference between parameter estimation which includes system identification and learning in machine learning perspective? Let say the model is y= Ax. x is the input and y is the output. In estimation, I have seen parameters to be estimated, maybe in this case it is A and the samples are also estimated (unsure) What is learning then?


title: Independence of effects in meta-analysis and meta-regression

text: In meta-analysis, I know that an important assumption is independence. Clearly, if we want to estimate a pooled effect size, its best if each effect comes from an independent sample of observers. Otherwise, the summary effect will be overly influenced by samples that contribute more than one effect size to the model. **However** , as an experimental psychologist, many of the comparisons I am interested in differ at the within subjects level (for instance, the same sample of subjects may have been tested on two sets of stimulus types). My question, therefore, is whether moderator-analyses can be conducted where k represents conditions, rather than samples. In other words: **If I had two effects for every sample, one for stimulus 1 and one for stimulus 2 -could I perform a moderator analysis on stimulus type, given that this varies at the within subjects level?** Note the following things- 1) I wouldn't use this "between condition" model to estimate the "overall" pooled effect size because of the issues I mention in the first paragraph. 2) I have seen published examples of such a "between condition" analysis, but the consequences/ limitations/issues involved are rarely discussed. Additionally, If anyone has any references/ examples as to why this is ok/ not ok I would be very grateful. I appreciate your time.


title: Getting expected value of future value with time varying data (credit card revolving and fee data) . Customer lifetime value

text: I have a credit card data and that contains monthly amount of revolving and amount of fee for each customer. As a bank perspective, I want to get the expected value of future revolving amount and fee amount in total (by using discount factor in infinite time horizon). I want to get the expected future value of revolving and fee respectively. The one way of my idea is that E(future) = Sum [ {prob(revolving at time t)*amount of revolving at time t} / {(1+discount)^t} ] Very simply we can think (revolving case) probability of revolving at time t = x/n (where x is number of revolving in data and n is total month in data) amount of revolving at time = mean of amount of revolving in data (constant) Actually I want to estimates the present value of expected future values. but above way is too naive. Is there any good idea that I can estimate the future value of revolving and fee values with monthly timeseries data ? (I have other variables like demographic, creditcard usage history and others). Thanks. I prefer R for program. Any suggestion with R package ?


title: The number of coin tosses needed if the proportion of heads is to lie within 0.05 of p with probability at least 0.9?

text: There's a question I'm not really sure if I did it right or even understand what its trying to say. > There is a coin which produces heads with an unknown probability $p$. How > many times should we throw this coin if the proportion of heads is to lie > within $0.05$ of $p$ with probability at least $0.9$? Hint: Answer is not > complete if it relies on $p$ and do not worry about continuity correction. This question comes from a chapter about binomial normal approximation. So far I know I have a $\mathcal B(n,p)$ and I need to find the number of $n$ so that $P(c/n \le 0.05*p) \ge 0.9$ with $c$ as the number of heads? I'm not sure if this is right. In turn $P( c \le 0.05np ) \ge 0.9$ So I'll normalize the binomial such that $$ P\bigg( Z \le \frac{0.05np - n*p}{np(1-p)} +0.5\bigg) \ge 0.9 $$ with $0.5$ as the continuity correction. Am I on the right track? I thought about using confidence intervals but we haven't discussed about that in the chapter or lecture.


title: Matroid techniques for identifying multicollinearity

text: Reading Belsey (Conditioning Diagnostics 1991) he mentions on page 35 what was then a new technique of using the flats of the Matroid of dependencies in X to identify multicollinearity see Greene below). http://www.sciencedirect.com/science/article/pii/016794739090108T He describes the technique as very promising but is too new at the time of publication for him to critique it meaningfully. He does note that that it represents an entirely new direction for research into identifying multicollinearity. However there seems to be little follow up literature. All I could find was the application used here: http://www.plosone.org/article/info%3Adoi%2F10.1371%2Fjournal.pone.0034410 Just wondered if anyone has experience using this technique and what it's pros and cons are relative to other multicollinearity measures such as the condition number? Also if the technique is valid why has it been cited so rarely in the literature?


title: Can I use Kolmogorov–Smirnov test on my Data?

text: I am not good in statistics so I desperately need your help. So I have this dataset of distributions, and I want to know if I can use the KS-test on it. the Idea is saying that the feature's distribution 1 in case 1 is different than its distribution 2 in case 2. For example, let s take the feature "size" in two cases: case1, case2 The distribution 1 (in case 1) looks like this:               [0,0,0,0,0,0,132,33,1200,0,0,98,208,56,0,0,0,....]      The distribution 2 (in case 2) looks like this:               [52215,2132,933,11200,0,0,13245,4208,309,0,34000,0,....]      and so on, each number represent the total size in one second, and the null hypothesis, is that **distribution 1 and distribution 2 follow identical distribution** so the point is rejecting it by having a less than 1% as a p-value (that s what I understood please correct me if I am wrong) I read that KS-test is applied on continuous distributions, is the one I have continuous?? how to know if your distribution is continuous? If I can't apply the KS, what else can I apply? mentioning that I work with Python..


title: Weighting more recent data in Random Forest model

text: I'm training a classification model with Random Forest to discriminate between 6 categories. My transactional data has approximately 60k+ observations and 35 variables. Here's an example of how it approximately looks like.                _________________________________________________     |user_id|acquisition_date|x_var_1|x_var_2| y_vay  |     |-------|----------------|-------|-------|--------|     |111    | 2013-04-01     | 12    | US    | group1 |     |222    | 2013-04-12     | 6     | PNG   | group1 |     |333    | 2013-05-05     | 30    | DE    | group2 |     |444    | 2013-05-10     | 78    | US    | group3 |     |555    | 2013-06-15     | 15    | BR    | group1 |     |666    | 2013-06-15     | 237   | FR    | group6 |      Once the the model is created, I'd like to score observations from the last few week. As there have been changes to the system, the more recent observations will resemble more closely the environment of the current observations that I'd like to predict. Hence, I want to create a weight variable so that the Random Forest would put more importance on the recent observations. Does anyone know if the randomForest package in R able to handle weights per observation? Also, can you please suggest what is a good method for creating the weight variable? For example, as my data is from 2013, I was thinking that I can take the month number from the date as weight. Does anyone see a problem with this method? Many thanks in advance!


title: Regression layer in convolutional neural network

text: I am trying to do regression with a deep convolutional network. I am using the code from http://deeplearning.net/tutorial/lenet.html Currently, the network is using logistic regression as its last layer which I want to remove. If I simply remove the layer               layer3 = LogisticRegression(input=layer2.output, n_in=500, n_out=10)      from the network, I will need to specify a cost instead of:               cost = layer3.negative_log_likelihood(y)      How can I specify an appropriate cost for doing regression? Each of my outputs should be 0 to 1, but not necessarily sum to 1.


title: Central Limit Theorem: Likelihood multiplication

text: Been watching this video by Tom Minka on Expectation Propagation (http://videolectures.net/mlss09uk_minka_ai/). At about 19:12, he says that the reason the moment matching technique works is that when you multiply a lot of likelihood terms together, the resulting posterior distribution is compact and its variance is very small and it will look more Gaussian (assuming IID). At some point in the video, he says this is because of the central limit theorem. However, However, I am unable to see how the CLT applies to multiplying such likelihoods and why this should approximate a compact Gaussian. Would appreciate any thoughts/suggestions on understanding this.


title: What is a "subject" in ICC?

text: I have been working in R and have finally got some code working regarding intraclass correlation coefficients. In the output I get something called "subjects". I thought it was to do with the number of columns of data but for a file with 20 columns I get 3 or 4 subjects depending on the data set. I have very little knowledge of R and I'm not great at statistics talk so would someone be able to explain in simple terms what a subject is when it is given in the output in R having done a single score intraclass correlation please? Thanks!


title: Political preference as independent variable - fixed or random effect?

text: I'm having a discussion with my colleague about an exam question. Let's say we have a questionnaire in which, among others, the preference for a political party is asked for. How would we model this, do we use a fixed or random effect? I think we should go for a random effect as it's a questionnaire and we expect the effect to change with the sample demographics (age, sex,...).


title: Popularity and Difficulty score calculations within groups and sub groups

text: I have an examination application that carries questions of various difficult groups: a) Basic b) Standard c) Advanced Based on data generated by users taking these questions I am trying to create a score of popularity and difficulty for the questions. I figured that three factors work towards this: 1) Number of times the question was taken (users can chose from a variety of topics) 2) Number of times the question was answered correctly (a question that gets moderately answered tends to be more 'popular') 3) Time spent on the question (indicator of difficulty AND relevance for the user, an assumption that the user wouldn't spend much time on a question they know nothing about) 4) Number of times the question was attempted (again indicates whether the question is actually relevant and creates an interest in the user) My Question: a) Can I evaluate these scores across the difficulty groups as a single score or do I need to evaluate it within the individual groups? b) Currently I calculate application averages of **attempts, time , correct attempts** and then divide each questions average time with the application average to calculate the score? What kind of averages or statistical formulas should I be considering for this kind of problem?


title: Understanding the role of the chi-squared distribution in the confidence interval for the variance

text: In my textbook they have this inequality: $$ \chi_{1-\frac{\alpha}{2}}^2 < \frac{(n-1)s^2}{\sigma^2} < \chi_{\frac{\alpha}{2}}^2$$ which later becomes this statement: $$\frac{(n-1)s^2}{\chi_{\frac{\alpha}{2}}^2 } < \sigma^2 < \frac{ (n-1)s^2}{ \chi_{1-\frac{\alpha}{2}}^2}$$ Now I know the whole idea is to find the confidence interval for $\sigma^2$ the variance, but I was wondering if the distribution for the variance is normal. I also don't understand why the chi square is squared. When I look at the picture in the book that shows a right skewed graph with the chi squares labeled (i.e. $\chi_{0.95}^2 = 4.575$ and $\chi_{0.05}^2 = 19.675$), I get the impression that I'm looking at something similar to $z$ scores. What are these chi squares? Do they represent the number of standard deviations away from the mean? I have also asked this question on math.stackexchange section but I haven't heard back from anyone.


title: What is the probability we have pairs of teammates that are teammates again in a second round?

text: There are 42 actors in our movie set. We have 2 rounds of random assignments to teams of 3. In each round, we will randomly divide the actors into 14 total groups of 3 persons each. (A) Suppose that we have that A and B are teammates in Round 1. What's the probability that they are teammates again in Round 2? (B) What's probability that no pair of teammates in Round 1 are teammates again in Round 2? I tried to use bayes formula for part A but it appears that the two events are independent of one another. So therefore, is it the case that the given information isnt needed? for b, I tried to do 1-P(at least one pair of teammates), does this approach work? Thanks! My work for A: P(teammates in round 2|teammates in round 1) = P(teammates in round 2 and teammates in round 1)/P(teammates in round round 1) = P(teammates in round 2) (since these are independent of each other). The probability of being teammates in round 2 is just $\frac{1}{14} * \frac{1}{14}*14$, where we have 14 ways of arrangements. So the final probability for part A is $\frac{1}{14}$


title: deletion residuals in linear

text: > **Possible Duplicate:**   >  Deletion residuals The $i^{th}$ deletion residual $e_{−i}$ is defined as $e_{−i}=y_{i}−X^{⊤}B_{−i}$ where $X^{⊤}$ is the $i^{th}$ row of the design matrix $X$ and $B_{−i}$ is a column vector of least square parameter estimates calculated without the $i^{th}$ observation. Write some annotated R code to calculate the deletion residuals when the linear model $y_{i}=B_{0}+B_{1}X_{i}+B_{2}X_{2i}+E_{i}$ is fitted to the data in the file quadratic.txt. By drawing suitable plot, comment on the distribution of these deletion residuals. should I write what the exactly distribution?


title: Active set methods for group lasso?

text: When the solution is extremely sparse, probably the fastest method for solving LASSO regression is least angle regression, which starts from an all-zero solution and adds nonzero elements to the parameter vector one at a time. I am trying to solve a group lasso problem: $ \min_{\beta} \frac{1}{2}||y-X\beta||+\lambda\sum_g||\beta_g||_2$ when the parameter space is extremely large, the columns of the design matrix $X$ are extremely correlated (basically a convolution matrix) and the optimal $\beta$ is known to be extremely sparse. Are there active set methods for the group lasso analogous to LARS? I know that the full regularization path is no longer piecewise linear like it is for lasso, but even a good approximate method would work well.


title: How to identify structural change using a Chow test on Eviews?

text: I have this little problem and I would appreciate some help. As part of my master thesis, I have to identify a trend in a univariate (GDP) time series for different countries. I have to separate the trend and the stochastic element in it for each country. I have managed to do so by doing: variable c @trend // for each country. And then running a AR(1) on the residuals // for each country. However, now I need to identify structural breaks in one of these countries. I've been reading and searching all over the internet and books and I've found that the test most people use to identify these structural changes is the Chow Test. I know how to run the test, but I have't been able to figure out how to interpret the results, and decide whether there is a structural break or not. Here there is an example of the results: ![enter image description here](http://i.stack.imgur.com/RL9Lz.jpg) What puzzles me the most is the fact that, regardless the point I choose to break the series, I always get Prob. F(2,47) 0.0016 //or any very significant value, with the same degrees of freedom. Can someone please help me understand how I should interpret these results in order to identify where the breaks lie?


title: Grading a range of data

text: I am new to data mining and I am exploring ways to divide a range of data into bins. It would be much easier to search for the algorithms if I knew the scientific term for this job (ranking maybe?). Here is what I desire: I have a data within a defined range, i.e. [min, max]. I know the first order statistics of the data like mean, variance, median, etc. I also plotted a histogram of it. My situation is somehow similar to assigning letter grades to students after an exam by examining the distribution of the overall exam grades. I want to have grades for the data from 1 to 5, with the first bin containing the smallest data values and 5 containing the highest values. I could do it uniformly by dividing my data range into 5 equal pieces but that won't be satisfactory. I would like to do this by inspecting my data distribution. I will appreciate if you could suggest me some known algorithms to do that.


title: How exactly are standard errors of coefficients (factor levels) computed in a two-way ANOVA?

text: In one-way ANOVA ($Y_{ij} = A_i + \epsilon_{ij}$) the standard errors of coefficients are computed pretty easily: $$SE_{\widehat{A_{i}}} = \sqrt{\text{mean residual sum of squares} \over \text{size of group i}}$$ How exactly are the standard errors of coefficients (factor levels) computed in two-way ANOVA? Is there similar expression? I did some research but no result. I've been trying to find out how is the variance-covariance matrix defined for the linear models (quite a basic thing I thought), but with no result! **I highly prefer direct expression** like the above with sums and indices over the matrix multiplication stuff. Thanks! :-) Here is example test code for two-way ANOVA               set.seed(2)     y = c(rnorm(200, 20, 4), rnorm(300, 30, 4), rnorm(400, 40, 4), rnorm(500, 50, 4))     cat1 = as.factor(c(rep(1, 200), rep(2, 300), rep(3, 400), rep(4, 500)))     rand_order = sample(1:length(cat1))     cat2 = cat1[rand_order]     cat2y = c(rep(1, 200), rep(-2, 300), rep(3, 400), rep(-4, 500))     y = y + cat2y[rand_order]     summary(lm(y ~ 0 + cat1 + cat2))


title: Repeated measures interaction betwen SNP and treatment with R

text: This is my first question is this forum. I have the a problem that I want to solve with the `lme` function of the `nlme` package of `R`. The goal of my study is analyse if the cholesterol level (quantitative) change in a series of patients depending of 3 different treatments and the genetic profile (SNPs, a factor variable with three categories). We have two repeated measures of the cholesterol (before the intervention and 1 year after). The general question is to know if some of the treatment "works" (change in the cholesterol) for some of the 3 categories of the SNP. My analysis in `R` is the next, but I don't know if it is correct:               w.lme.1 <- lme( colesterol ~  tiempo * treatment + SNP ,                   random = ~1 | ID , data=xx, na.action=na.omit, method="ML")          w.lme.2 <- lme( colesterol ~  tiempo * treatment * SNP ,                   random = ~1 | ID , data=xx, na.action=na.omit, method="ML")          anova(w.lme.1, w.lme.2)      where:   * ID , identification of the individuals   * time, 2 factors : baseline time and 1-year after   * treatment, 3 levels   * SNP, 3 levels R output:               > head(xx)        ID colesterol tiempo treatment SNP     1   1     312.47      0         1   1     2   1     221.61     12         1   1     3  10     221.33      0         1   1     4  10     227.67     12         1   1     5 100     229.07      0         2   1     6 100         NA     12         2   1     > str(xx)     'data.frame':   594 obs. of  5 variables:      $ ID        : int  1 1 10 10 100 100 101 101 102 102 ...          $ colesterol: num  312 222 221 228 229 ...      $ tiempo    : Factor w/ 2 levels "0","12": 1 2 1 2 1 2 1 2 1 2 ...          $ treatment : Factor w/ 3 levels "1","2","3": 1 1 1 1 2 2 1 1 3 3 ...      $ SNP       : Factor w/ 3 levels "0","1","2": 2 2 2 2 2 2 1 1 1 1 ...


title: Epoch vs Incremental Update over Entire Dataset

text: In the context of neural networks, what is the difference between training like this:               foreach observation in dataset:         compute_error(observation)         update_weights()      and               foreach observation in dataset:         compute_error(observation)     update_weights()      Both train over the entire dataset.


title: Prediction with GLS

text: Let's say I build a Generalized Least Squares model. I follow the standard procedure and first estimate a LM model. Then I create an error-response covariance matrix based on the residuals of this model. Now I build an LM model again only this time I specify weights based on the error-response covariance matrix. Now suppose I want to predict with the GLS model out-of-sample to test for model stability. I want to confirm that I can simply perform a prediction using the coefficients estimated by GLS and there is no need to furnish weights anymore (especially since in a prediction scenario where residuals are not available the error-response covariance matrix cannot be generated). Follow-up question: We proceed to score on test data with coefficients from the training data. (The dimension of the test data consists of a cross-section of N individuals and T observations.) We would like to produce consistent standard errors. Therefore, instead of calculating standard error of the estimate in the OLS fashion, we weight the residuals (see below) by a "GLS weights" vector: OLS calc of SEE: `sqrt( sum( ( residuals from linear model ) ^ 2 ) ) / residualDegreeFreedom )` GLS calc of SEE: `sqrt( sum( ( residuals from linear model) ^ 2 * glsWeight ) ) / sum( glsWeight ) * length( glsWeight ) / residualDegreeFreedom )` "gls weight" is a vector calculated in the usual way as the inverse of the variance of the residuals of each cross-section at a date (i.e. a vector of length T). However, here I am using the residuals from the _test_ data as opposed to the training data (indeed this is required otherwise the dimension of out-of-time residuals would not match the dimension of GLS weights vector). What is counter-intuitive is that if I want to measure the SEE of the GLS model out-of-sample on _one_ individual, I am required to score _all_ individuals out of sample (otherwise constructing the GLS weights vector would be impossible since there is no variance of residuals). Question is - Am I required to use the GLS weights when calculating the SEE out-of-sample, or can I simply use the OLS calculation of SEE?


title: (Stochastic) optimization with multivariate predictive distribution

text: Suppose I have a multivariate predictive distribution (dimension N). The distribution cannot be described in closed form with a mu vector and associated covariance matrix. However, I can sample from the distribution various scenarios with associated probabilities. One can imagine sampling 100,000 scenarios to construct a matrix where each row represents the joint distribution for a scenario, and a probability is associated with each scenario. Now let's say I have a weight vector and I am performing some optimization involving the distribution (for example, a mean-variance optimization) with some constraints. Naturally, for any given draw from the distribution we can search for the optimal weights that maximize the mean-variance function: weights * sampled return - weights(transposed) * E * weights. However, the optimal weights for one draw will in general not correspond to the optimal weights from another draw. Question - what is the optimal choice for the weight vector? I believe the answer may involve a stochastic optimization procedure, however, I am not familiar with this space and would like to validate this. There are also natural connections to the literature on robust optimization but the solutions here require closed form expressions for the multivariate distribution. I have explored re-sampling which has an analogy to boot-strapping. This involves constructing optimal weights for many samples from the predictive distribution and somehow averaging or combining the weights. However, this approach is flawed and lacks a theoretical basis. How might I proceed to tackle this problem with stochastic optimization or some other technique?


title: Quantifying the loss of efficiency incurred by a pseudo-likelihood estimator relative to the true maximum likelihood estimator

text: I'm working on a couple of problems where I use pseudo-likelihood (a.k.a. composite likelihood) as an estimation method since full maximum likelihood is not feasible. To give some background without getting into very specific details, in both projects the data are iid samples of clustered binary data such that there are dependencies within each cluster. That is, I have iid observations of $Y_{1}, ..., Y_{n}$ such that $Y_{k} = \\{ Y_{k1}, ..., Y_{km_{k}} \\}$ where $m_{k}$ can be rather large (say, 100). Since there is no simple a prior choice for the distribution of multivariate binary data, I attack the problem by modeling the bivariate distributions and producing a pseudo-likelihood by summing over the pairwise likelihoods. This is sensible because the likelihood of a bivariate binary vector is a simple function of the cell probabilities and can be conveniently re-parameterized by the marginal probabilities and the odds ratio, giving interpretable parameters in terms of the mean and dependence structure. Specifically, I maximize $$ {\rm CL}(\Theta) = \sum_{k=1}^{n} w_{k} \sum_{(i,j)}\log \Big( P(Y_{ki}, Y_{kj} | \Theta) \Big ) $$ as a function of the parameter $\Theta$. The index $(i,j)$ is meant to represent all pairs of observations within cluster $k$. The weights $w_{k} = m_{k}/\binom{m_{k}}{2}$ are incorporated so that a particular cluster doesn't have an artificially large influence due to the pairwise expansion. It's a fact that pseudo-likelihood estimators are, at best, as efficient as full maximum likelihood estimates, and are usually less efficient (clearly). In other work I've read whenever pseudo-likelihood methods are used there is some discussion/results attempting to quantify the loss of efficiency. This leads me to my question: Since I can't feasibly calculate the maximum likelihood estimate in terms of my model parameters (except in simple cases that are actually sub-models captured exactly by the composite likelihood estimator, such as independence), how exactly can I quantify this loss of efficiency? Would it be at all convincing to anyone if I were to present estimation results when the true model is and is not a subset of the space covered by the composite likelihood estimator (since it is simple enough to simulate data with higher order dependencies than those captured by the pairwise likelihoods) and show that not too much precision is lost? My initial thought was to calculate the asymptotic relative efficiency of the maximizer of ${\rm CL}(\Theta)$. I have an analytic expression for the asymptotic covariance matrix of the pseudo-likelihood estimator which I can numerically compute if I must (but that will require me to calculate the hessian, ugh). But, I can't make much headway with the analogous quantity for the MLE, so that's dead in the water, unless I'm missing something. Thanks for any thoughts.


title: What will be the effect on the d.f of $\chi^2$ distribution

text: Cochran's Theorem says that > The conditional distribution of $\chi^2=\sum_{i=1}^nX_i^2$ under $m(<n)$ > independent homogeneous linear constraints: > $$\begin{align}a_{11}X_1+\cdots\dots +a_{1n}X_n&=0\\\\\vdots > \\\a_{m1}X_1+\cdots\dots +a_{mn}X_n&=0\end{align}$$ is $\chi ^2_{n-m}$, > provided $X_1,\dots,X_n$ are _i.i.d_ $N(0,1)$ variables. It is clear that if the constraint is linear then the _dergees of freedom_ (d.f) diminise.If we consider non-linear (like quadratic or cubic) constraint what will be the effect on the d.f of $\chi^2$ distribution.


title: Negative coefficient in ordered logistic regression

text: Suppose we have the ordinal response $y:\\{\text{Bad, Neutral, Good}\\} \rightarrow \\{1,2,3\\}$ and a set of variables $X:=[x_1,x_2,x_3]$ that we think will explain $y$. We then do an ordered logistic regression of $X$ (design matrix) on $y$ (response). Suppose the estimated coefficient of $x_1$, call it $\hat{\beta}_1$, in the ordered logistic regression is $-0.5$. How do I interpret the odds ratio (OR) of $e^{-0.5} = 0.607$? Do I say "for a 1 unit increase in $x_1$, ceteris paribus, the odds of observing $\text{Good}$ are $0.607$ times the odds of observing $\text{Bad}\cup \text{Neutral}$, and for the same change in $x_1$, the odds of observing $\text{Neutral} \cup \text{Good}$ are $0.607$ times the odds of observing $\text{Bad}$"? I can't find any examples of negative coefficient interpretation in my textbook or Google.


title: How to modify RankBoost to maximize area under recall-precision curve instead of AUC?

text: Using the WeakLearn algorithm from the original RankBoost paper, how do you set the optimal threshold to maximize AU-RPC (instead of AUC)? And, once that threshold is set, how do you calculate the AU-RPC for use in calculating alpha? This paper suggests that the "area under the recall-precision curve [is] a much better metric [than AUC] when working with highly skewed data sets", which is my current application. The paper proposes a modification of RankBoost.B that purports to maximize AU-RPC. However, the authors are using learners based on first-order clauses, while I need to use the originally proposed WeakLearn learners. How can I apply the same idea with the original learner algorithm?


title: Frequentist statistics references for someone well versed in modern probability theory

text: Coming from a rigorous background in analysis and modern probability theory, I find Bayesian statistics straightforward and easy to understand, and frequentist statistics incredibly confusing and unintuitive. It seems that frequentists are really doing bayesian statistics, except with "secret priors" that aren't well motivated or carefully defined. On the other hand, a lot of great statisticians who understand both perspectives ascribe to the frequentist perspective, so there must be something there that I just don't understand. Rather than giving up and declaring myself a Bayesian, I'd like to learn more about the frequentist perspective to try to really "grok" it. What are some good references for learning frequentist statistics from a rigorous perspective? Ideally I'm looking for definition-theorem-proof type books, or perhaps hard problem sets that, by solving them, I would gain the right mindset. I've read a lot of the more "philosophical stuff" one might find searching the internet - wiki pages, random pdfs from .edu/~randomprof sites, etc - and it hasn't helped.


title: How do you choose the timeframe for a prediction to be tested?

text: To (hopefully) clarify the question a bit more if you have a broad working hypothesis that x influences y, on a causal/predictive basis (e.g. if x goes up in time t, then you would expect y to go up in time t+1, assuming the working hypothesis is true). In traditional time-series analysis the length of time that the period t represents for t and t+1 tends to be equal. Can alterations in the periods that each represents, be easily altered for the purposes of analysis? Also if we take the example of t to be a day, is it possible to distinguish differences in predictivity for various time frames, e.g. x goes up today, and over the course of the next 3 days y goes up, or over the course of the next 5 days y goes down (what happens if they are both significant?). How does one go about detecting the time-frame of prediction appart from brute-forcing every possible combination of before and after timeframes (e.g. testing if x has gone up over the course of the last m days what happens to y over the next 1 to n days).


title: Correlation and non-normal distributions

text: I have to work with bimodal and normal distributions. Starting from a set of given bimodal distribution I have to draw a set of random vectors. Some of these random vectors must have given correlation coefficients. I thought to draw a vector $x_1$ of random values from one of the bimodal distributions, and then mix with some gaussian noise $\epsilon$, in order to get a new vector $x_2$ correlated with the first. But I need to know a priori the correlation coefficient. Moreover I think I can't use pearson's correlation since mean and variance doesn't make sense for a bimodal distribution. Any advices?


title: Finding multivariate clusters with survey data (in R)

text: I'd like to conduct a multivariate cluster analysis on data from the American Community Survey's PUMS microsample (individual level records). I've only performed cluster analysis before when there are no errors (or imagined to be no errors) in the data. Has anyone tried to combine Thomas Lumley's excellent survey analysis package for R with any of the cluster analysis packages? `Pvclust`, for instance, provides p-values to identify statistically significant clusters. However, I do not know how to feed the standard errors from the survey package into a search for multivariate clusters. The problem is compounded by the presence of categorical variables relevant to my research.


title: implementing Lasso with BIC

text: Do you have an R code to implement `Lasso` with `BIC`? Note that there is an R package called `parcor` to implement **lasso** and **adaptive lasso with Cross validation**. But what I am looking for is the **lasso with BIC**. Can you please help me?


title: How does "ward" clustering (in R's hclust function) work?

text: A simple example: `plot(hclust(dist(c(1:3)),method = "ward"))` I would like to know which calculations (in R) can reproduce the distance of 3 from {1,2} to be 1.67 ![enter image description here](http://i.stack.imgur.com/XXMs2.png) Thanks.


title: Modelling worst, average and best case capacity of a system

text: A system has an incoming request rate of N requests per second. The system hosts a pool of x workers where x >= N. Each worker can complete a single request in approximately t seconds. The maximum parallel requests which can be processed by the system is x. But sometimes a worker can get delayed indefinitely and when it exceeds T seconds the current request is aborted and a new request processed _(t < T)_. This could happen together for 1% - 100% of the workers for several seconds. During those T seconds the worker is blocked and cannot process any other request. These unprocessed requests are dropped and never processed again. To avoid being indefinitely blocked, the system aborts the current request to start process a new request. On an average n requests are aborted every hour. How do you model this information to find out the following ?   1. The number of dropped requests as a function of requests aborted.   2. The effect of tuning timeout(T), pool size of workers(x), average process time(t) on reducing or increasing the number of dropped requests or it's inverse, the number of processed requests ?   3. The worst, average or best case capacity of the system.


title: Adjust a large set of Spearman correlation analyses for multiple testing

text: > **Possible Duplicate:**   >  Threshold for correlation coefficient to indicate statistical significance > of a correlation in a correlation matrix ### Context I am doing an exploratory study to investigate the relationship between a drug (actually measured in two ways - by direct and indirect methods) and 15 various parameters. There are three groups of individuals, on which I compute correlation statistics separately ($n=11$, $n=11$ and $n=24$ - the first two are paired). Obviously, the number of statistical tests is large, and I would like to control for multiple testing. I am using Spearman correlation as data are limited and do not (always) follow a normal distribution. I am working in R. ### Questions   * How should I control for multiplicity of testing?  This is an hypothesis generating study, and I would rather prefer to use less conservative methods. I have been looking into Benjamini–Hochberg's procedure (known under R as `"BH"` or `"fdr"`) which allows to control the false discovery rate. However, I am not sure if I would violate dependence assumptions. And I am not sure if this could be used in correlation statistics. Perhaps I should not adjust for multiple testing at all. In my mind, if both direct and indirect methods give the same association, false positives are highly unlikely.


title: Probability of an event that is not measureable

text: We know from measure theory that there are events that cannot be measured, ie they are not Lebesgue measureable. What do we call an event with a probability that the probability measure is not defined on? What types of statements would we make about such an event?


title: Calculating precision-recall curve from error matrix in MATLAB

text: How is a precision-recall curve calculated from an error matrix and vector of labels? The error matrix is of size NxN, and contains distances between the features associated with the respective row and column (the diagonal is consequently 0). The labels are in a 1xN vector, describing the class corresponding to the rows/columns of the error matrix, e.g. [1 1 1 1 2 2 2 ... 20 20 20 21 21 ... C-1 C-1 ... C C C], where C is the number of classes. I've used the two equations below on all rank results for each class to end up with precision and recall values for each class.               Precision = TP / (TP + FP)     Recall = TP / (TP + FN)      My question is (if what I'm doing is correct): is how are all the class precision-recall values combined/averaged to obtain a single plot? Do I just average all precision values at a given value of recall? Or is it acceptable to average all precision values for a given rank, and average all recall values for a given rank? Thanks very much for any insight. The precision recall curves I get are shown in this image. ![enter image description here](http://i.stack.imgur.com/JhVkp.png) Code following (but probably unnecessary).               % returns the precision-recall points for plotting.      % 'dismat' is an NxN distance matrix containing the distance between each feature and every other feature     % 'labels' is an N-length vector that that corresponds to each row and column of 'dismat'     function [prec,recall] = calcPrecisionRecall(dismat, labels)          if size(labels, 1) ~= 1         labels = labels';     end          numClasses = length(unique(labels));          [sortedDistMat, idx] = sort(dismat);   % every column is sorted     recalledClasses = labels(idx);          recalledMatches = bsxfun(@eq, recalledClasses(2:end, :), labels);          classMatches = zeros(size(recalledMatches, 1), numClasses);     instancesPerClass = zeros(1, numClasses);     for ii=1:numClasses         classMatches(:, ii) = sum(recalledMatches(:, labels == ii), 2);         instancesPerClass(ii) = sum(labels == ii);     end          % we'll just increase the rank to get the points     prec = zeros(size(classMatches));     recall = zeros(size(classMatches));     sumClassMatches = cumsum(classMatches);     rank = 1:size(classMatches, 1)     for ci =1:numClasses         prec(rank, ci) = sumClassMatches(rank, ci) ./ (rank' * instancesPerClass(ci));         recall(rank, ci) = sumClassMatches(rank, ci) / ((instancesPerClass(ci) - 1) * instancesPerClass(ci));     end          % average all precision and all recall values (datapoints will therefore be avg prec--avg recall for a certain rank)     avgPrec = mean(prec, 2);     avgRecall = mean(recall, 2);      I previously asked this question on StackExchange (http://stackoverflow.com/questions/19824979/matlab-precision-recall-from- error-matrix) but was told to bring my question here.


title: Likelihood ratio test statistic

text: Let $X_{1},..,X_{n}$ be a random sample from a normal distribution with mean ${\theta}$ and variance 1. Let $H_{0}:{\theta}=0$ and $H_{1}:{\theta}{\neq}0$. Using likelihood ratios, show the critical region is of the form $exp(\frac{-n{\bar{x}}^2}{2}){\geq}k$ for some non-negative $k$. My work: taking the supremum of the likelihood function over non-zero ${\theta}$ gives $1/({\sqrt{2{\pi}}})^n$. So the test has critical region of the form $exp(0.5(x_{1}^2+...+x_{n}^2)){\geq}A$ for some non-negative $A$. How can I extract the square of the mean from the sum of squares? Thanks


title: How to test for significance if the variances are different?

text: If got data (value) of two types (cond true or false):                 > str(test.df)       'data.frame': 3208 obs. of  2 variables:        $ cond : Factor w/ 2 levels "FALSE","TRUE": 1 1 1 1 1 1 1 1 1 1 ...            $ value: num  31 25 21 29 18 41 15 7 33 6 ...      The mean and the variance (and also the number of cases) are different:                   > stats <- ddply(test.df, .(cond), summarize, mean=mean(value), var=var(value), n=length(value))         > stats            cond     mean      var    n         1 FALSE 17.33918 141.3199 3137         2  TRUE 25.91549 177.4499   71      So different variances means I can't use wilcoxon-test, right? If I plot the densities I see the following: ![densities](http://i.stack.imgur.com/s4Wt8.png) The Shapiro-test says they are not normal-distributed.                   > shapiro.test(test.df[test.df$cond==TRUE, 'value'])                Shapiro-Wilk normality test              data:  test.df[test.df$cond == TRUE, "value"]         W = 0.9583, p-value = 0.01901              > shapiro.test(test.df[test.df$cond==FALSE, 'value'])                Shapiro-Wilk normality test              data:  test.df[test.df$cond == FALSE, "value"]         W = 0.9348, p-value < 2.2e-16      No normal distribtion means I can't use Welch-test, right? So how can I test if the difference of the means is significant?


title: How do you solve a Poisson process problem

text: During an article revision the authors found, in average, 1.6 errors by page. Assuming the errors happen randomly following a Poisson process, what is the probability of finding 5 errors in 3 consecutive pages? Please explain your methodology, as the main purpose of the question is not getting the answer but the "how to". Please consider this as homework if you will. Thanks!


title: Integrating p-values (or z-scores) with missing values

text: I would like to apply a Stouffer (or Fisher) integration on a list of p-values or z-scores. However, some of them contain missing values. E.g. (R snippet):               a<-c(0.05,NA,0.04)     b<-c(0.03,0.13,0.01)      What is in your opinion the best approach, in order to make the integrations of a and b comparable? 1) simply remove NA data and calculate. Both Fisher and Stouffer methods take into consideration the number of values to integrate 2) consider the NA data as "not significant" (sounds like a bad solution) Any suggestion is much appreciated. Thanks!


title: Do I have to cite 'lme4' when using 'lmerTest'?

text: I am using `lmerTest` to run linear mixed models (LMM) to get the p-values. However, in the articles I have written, they will show both `lme4` and `lmerTest`. Then I am not sure if I also have cite `lme4` when I use only `lmerTest` to run LMM models. Do you have any idea about this?


title: Importance of Time Features

text: if you have a time series and you want to do some predictions, what **time feature** should you use ? lets say we are trying to predict _how many people visit a certain website_ , we have data for the visits for the last 2 years, what should we include as a time feature ? and if we use more features might it make our model less accurate ? **the features i can think of are the following:**               Date index(1-700)     Week number (1-53)     Day of the week (1-7)     Month(1-12)     Day of the month(1-31)     Day of the year(1-365)     Year


title: Moderator Sample Size Discrepancy in Meta Analysis

text: I'm coding a meta analysis and have come across a paper that seems to have significant variation in the number of participants who answered specific demographic questions. The total sample size is $N = 81$, but the number of participants used to calculate and report age ($n = 72$), years of career experience ($n = 74$), and gender ($n = 80$) varies. Is there a rule of thumb for whether/how to code this data for moderator analyses? I've checked the Cochrane Collaboration but have been unable to find anything. My concern is that by recording the sample size as 81 (the $N$ used for the effect size), that it will in some way skew my moderator analyses because for age, the actual n is much lower. Should I simply list the data as "missing" for this particular study? Can I transform it, or even use it as is?


title: Wrong parameters for GEV

text: I am doing some data analysis involving fitting datasets to a GEV distribution, but I'm getting some weird results. I'm using scipy, which uses MLE for fitting the parameters. My data is               [1.47, 0.02, 0.3, 0.01, 0.01, 0.02, 0.02, 0.12, 0.38, 0.02, 0.15, 0.01, 0.3, 0.24, 0.01, 0.05, 0.01, 0.0, 0.06, 0.01, 0.01, 0.0, 0.05, 0.0, 0.09, 0.03, 0.22, 0.0, 0.1, 0.0]      Using `scipy.stats.genextreme.fit`, I get the following parameters:               shape = -9.2409     location = 0.003991     scale = 0.03688      The Kolmogorov-Smirnov test (`scipy.stats.kstest`) gives               D = 0.4171     p-value = 2.1576e-08      And a plot of the histogram and pdf looks like this: ![histogram & pdf](http://i.stack.imgur.com/ku8SP.png) Now the problems:   1. The data are hydrological flows. If I use the GEV model to extract, e.g. the 50 year recurrence flow I get a ridiculously high value, something in the order of 1e17, so intuitively the model is incorrect.   2. The p-value is very low, indicating that the model parameters are incorrect, or that the model itself is inappropriate. However, I don't know of any other that is suitable for extreme values.   3. We have a black-box spreadsheet developed many years ago which calculates the GEV using method of moments, and it appears to give reasonable results for specific recurrence intervals. Unfortunately the person who developed it is long gone and I can't access the calculations to see what's going on. However, it implies that method of moments may give a better fit than MLE for this data. Am I using the wrong methods for this type of data, or is there something wrong with the way I'm applying them?


title: factorial with missing cells R

text: I'm trying to analyse data from a factorial experiment with missing cells. I've read this reference http://users.monash.edu.au/~murray/BDAR/ which recomends to treat the experiment with an one-way anova, and make contrasts to evaluate the main effects and interactions. following the book examples I did               contrasts(fdata$ALL) <- cbind(c(rep(1,46),rep(-1,45), rep(0,32)),c(rep(-1,46),rep(0,45), rep(1,32)))     AnovaM(aov(resp ~ ALL, fdata),     split = list(ALL = list("PIE" = 1:2, " PIE1 vs PIE2"=1, " PIE1 vs PIE3"=2)))      I'm having trouble generating contrasts in R, and knowing how to generate the correct contrasts to each interaction. I have six factors, five with 3 levels, and one with two, does any one know a systematic way to treat this kind of problem?


title: Goodness of fit tests robust against non-normality and estimated parameters

text: I'm struggling to find a goodness of fit test of the above. The non-parametric tests I have looked at (KS) seem to be unable to deal with estimated parameters - can someone help?


title: linear regression after cluster analysis (How much information did I lose)

text: I had dataset of the type $y;x_1, x_2, x_3, ... $ and I used some of my variables to perform a cluster analysis. Let's say one of those variables that I used was $x_i$ Then. After I had my clusters. I decided to run a simple linear regression of the type $y \hat{} x_1 + x_2 + ...$ ${\text{within each of the clusters}}$ ; and not surprisingly I noticed that variable $x_i$ is almost a constant. This is because the original variable had a range of $[0;10]$. But the new cluster has a range of $[4;4.5]$ (that's just an example). So my guess is that my variable lost explanation power after I performed my cluster analysis. I would like to hear opinions as to how to measure how much information power I lost and whether I should use variables at all after performing a cluster analysis with them.


title: Normality Testing of Productivity Data

text: I have a series of measurements from a construction production study I've been carrying out. I essentially measured the output rate of workers in terms of output/hr. The measurements varied each day and I obtained a series of daily averages. So x variable would be date and y-variable would be the output rate. I carried this study on over approximately a month, yet I'm required to carry out some statistical analysis on this data. I've tested for normality using a Q-Q plot which shows a linear result, so that I can start to calculate a confidence interval. My question lies into whether what I'm doing is valid or not? I'm quite rusty on stats and am at a bit of a loss as to what statistical analysis I should be carrying out on my data set. Any help or direction as to what sort of statistical analysis would be valid on my data would be great.


title: How to calculate the confidence interval of a regression prediction given a particular value for a binary predictor?

text: I have a regression predicting wage (dependent variable) from whether the participant is a college graduate (dummy variable, independent). I have the regression coefficients and their standard errors, the R-squared of the regression, and the covariance between the intercept and the coefficient on the dummy variable. How do I calculate a confidence interval for the average wage of a college graduate (given that dummy = 1)?


title: Find sudden change or high variation in a dynamic real-time stream of real continuous number

text: I have a stream of real continuous numbers (Figure. 1). The stream flows from **left** to **right**. At $t_1$ the first number entered the stream. Every time a new number enters the stream I check if there is a sudden change or high variation. How can I compute the `sudden change` or `high variation` without using any **threshold**. The following picture is **just** and example. ![enter image description here](http://i.stack.imgur.com/xMCbs.png) * * * What I was doing is simply compute a quantity $Q_x = |t_4 - t_5|$ and if $Q_x > \delta$ then there is a sudden change and a high variation. Using a threshold is wrong in my situation because there is no bound for the real numbers in the stream. Another idea came to my mind which is the following: $if\;|t_5-t_4| >|t_4-t_3| + \delta $ where $\delta = [0,1]$ then there is a sudden change. What do you think? Thank you.


title: Why do we have initial wiggle in the prediction of ARIMA model

text: I have seen many prediction of time series using ARIMA model in R forecast package, and the result always have some initial wiggle before reaching the stable level. What causes the initial wiggle? I can't attach a image because I need at least 10 reputation to post images. Basically the predicted value for the first few steps forecast is different from the long stable level, but it converges to it after a few steps.               install.packages("forecast")     library(forecast)     install.packages("Ecdat")     data(Tbrate,package="Ecdat")     attach(as.list(Tbrate))     auto.arima(pi,max.P=0,max.Q=0,ic="bic")     fit = arima(pi,order=c(1,1,1))     forecasts = predict(fit,36)     plot(pi,xlim=c(1980,2006),ylim=c(-7,12))     lines(seq(from=1997,by=.25,length=36), forecasts$pred,col="red")     lines(seq(from=1997,by=.25,length=36), forecasts$pred + 1.96*forecasts$se, col="blue")     lines(seq(from=1997,by=.25,length=36), forecasts$pred - 1.96*forecasts$se, col="blue")


title: Is there a good browser/viewer to see an R dataset (.rda file)

text: I want to browse a .rda file (R dataset). I know about the `View(datasetname)` command. The default R.app that comes for Mac does not have a very good browser for data (it opens a window in X11). I like the RStudio data browser that opens with the `View` command. However, it shows only 1000 rows and omits the remaining. Is there a good browser that will show all rows in the data set and that you like/use. Thanks.


title: Proper likelihood function in acceptance probability of Gibbs Sampler

text: I have a question about the acceptance ratio used when implementing a random walk M-H in a gibbs sampler to generate sample paths of an unobservable process. When computing the likelihood of a set of parameters, does the likelihood also include the likelihood of each jump in the latent process given $ \theta $ ? The following is a simplified motivating example. Data is generated from the following process: $$ c_{i,t} = \alpha + \beta y_{t} + \epsilon_{i,t} $$ where $$ dy_{t} = -\kappa y_{t} + N(0,1) $$ Thus, $ Y $, which is unobservable and independent of the constant term, follows a Ornstein–Uhlenbeck process. I'm trying to estimate $ \theta = (\alpha, \beta, \kappa) $ using the EM algorithm where samples of the OU process are drawn from a gibbs sampler using random walk M-H. When drawing samples, I start with $ Y^{(0)}_{t} = 0 $ for $ t = 0,1,2...T$. I'm going to draw a new candidate $y_t$ from $ N(Y^{(n-1)}_t, 2) $ with acceptance probability $$min (\dfrac{L(\theta | C, Y^{(n-1)}_{(-t)},Y_t = y_t)}{L(\theta | C, Y^{(n-1)}_{t})}, 1) $$ My question is: Are the likelihood functions in the acceptance probability a joint probability of observing all values of $ C $ and the probability of going from $ Y_{t-1} -> y_t -> Y_{t+1} $ or just of observing all values of $C$? In other words, does the likelihood function look more like (when restricting the likelihood function to values of $c_{i,t}$ for one time period) $$ \prod_{j=1}^N \phi(c_{j,t} - \alpha - \beta y_t) $$ or $$ \lbrace \prod_{j=1}^N \phi(c_{j,t} - \alpha - \beta y_t) \rbrace \cdot p(y_t | Y^{(n-1)}_{t-1}) \cdot p(Y^{(n-1)}_{t+1} | y_t)$$ I have a similar question for the likelihood function I'm trying to maximize in the M-step of the EM algorithm. I assume when finding the $\kappa$ that maximizes the likelihood function across all sample paths, I'm supposed to include the joint likelihood of jumping from $y_{t-1}$ to $y_t$ for each jump in each sample path?


title: What is the power of the regression F test?

text: The classical F-test for subsets of variables in multilinear regression has the form $$ F = \frac{(\mbox{SSE}(R) - \mbox{SSE}(B))/(df_R - df_B)}{\mbox{SSE}(B)/df_B}, $$ where $\mbox{SSE}(R)$ is the sum of squared errors under the 'reduced' model, which nests inside the 'big' model $B$, and $df$ are the degrees of freedom of the two models. Under the null hypothesis that the extra variables in the 'big' model have no linear explanatory power, the statistic is distributed as an F with $df_R - df_B$ and $df_B$ degrees of freedom. What is the distribution, however, under the alternative? I assume it is a non-central F (I hope not doubly non-central), but I cannot find any reference on what exactly the non-centrality parameter is. I am going to guess it depends on the true regression coefficients $\beta$, and probably on the design matrix $X$, but beyond that I am not so sure.


title: How to do 1-step Chow test?

text: Has anyone heard of 1-step Chow test? I've searched the keyword "1-step chow test" in Google but can not find any hits. What's the formula of this test?


title: Do orthogonal projectors change the statistics of i.i.d. generated data?

text: I am wondering the orthogonal projectors such as Hadamard transform change the statistics of i.i.d. generated data. Consider, i.i.d. generated data $\mathbf{X}$ of length $N$, $$\mathbf{X}_P=\mathbf{P}\mathbf{X}$$ where $\mathbf{P}$ is an orthogonal projector of size $N \times N$ and $\mathbf{X}_P$ is the projected data of length $N$. It is explicit that $\mathbf{P}$ doesn't change the covariance matrix of i.i.d. data, but I am doubtful whether the projected data is still i.i.d. and the shape of distribution is the same, e.g., if $\mathbf{X}$ i.i.d. is generated from Laplacian, then $\mathbf{X}_P$ is still i.i.d. Laplacian? Thanks a lot in advance


title: Aggregating results from Arima runs R

text: /edit: To clarify: The mtable function from the memisc package does exactly what I need, but unfortunately does not work with arima models. Similar to this question: I have multiple Arima models, some of which I've also fit with dependent variables. I'd like an easy way to make a table/graph of the coefficients in each model, as well as summary statistics about each model. Here is some example code:               sim <- arima.sim(list(order = c(1,1,0), ar = 0.7), n = 200)          ar1<-arima(sim,order=c(1,1,0))     ar2<-arima(sim,order=c(2,1,0))     ar3<-arima(sim,order=c(3,1,0))     ar4<-arima(sim,order=c(2,2,1))          #Try mtable     library(memisc)     mtable("Model 1"=ar1,"Model 2"=ar2,"Model 3"=ar3,"Model 4"=ar4)     #>Error in UseMethod("getSummary") :      #  no applicable method for 'getSummary' applied to an object of class "Arima"          #Try  apsrtable     library(apsrtable)     apsrtable("Model 1"=ar1,"Model 2"=ar2,"Model 3"=ar3,"Model 4"=ar4)     #>Error in est/x$se : non-numeric argument to binary operator


title: Simplest test for consistency of a rating scale

text: I want to explain to a non-novice doctor the best way to validate consistency of a measure across a likert-scale. The variable describes an average dose across 5 levels of patients (lower to higher). Whats the easiest way to validate that the dose statistically increases from level 1 to level 2 to level3 ..to level 5?


title: Normalization with error in denominator

text: I am trying to come up with a proper way to normalize my data. Using a microscope I want to count the percentage of green cells in a population. However, only ~0.1% of all cells are green. I decided I want to count 30 green cells, this would mean that I need to count roughly 30,000 cells to find 30 green ones. However, I do not want to manually count 30,000 cells. So I came up with a short cut. I know the exact area of a single frame under the microscope. So if I know how many cells are in a single frame, and I keep track of how many frames I imaged, I do not need to count individual non-green cells. To estimate the number of cells/frame I counted 10 random frames, this allows me to estimate an average (+/- stdev) number of cells per frame and this gives me the denominator I need.   However, I am stuck at the next step. So I count green cells and have an extrapolated total number of cells. When reporting the % of green cells, I would normally # green cells/total cells. But now, I need to keep the error of the estimate into account.   I am not sure how I would do that. Furthermore, I am treating these cells with different drugs to see whether the number of green cells changes. What statistics should I use, or does the variation of the denominator not matter?   thanks!!


title: Non-algebric curve-fitting along weighted pointcloud (if possible using python)

text: I have a list of weighted 2D points taken from symmetry analysis of a human back surface. I am supposed to find the "midline" representing the most likely path describing vertebrae location (actually, spinal process of the vertebrae). EDIT: A representative dataset is as follows:               [[ -0.7898176   -3.35201728   4.36142086]      [  2.99221402  -3.35201728   1.11907575]      [  6.97475149  -3.35201728   2.4320322 ]      [ -4.82443609  -2.35201728   0.6479064 ]      [ -1.32418909  -2.35201728   1.88004944]      [  0.07067882  -2.35201728   1.10982834]      [  3.09169448  -2.35201728   1.8557436 ]      [  7.10399403  -2.35201728   2.03906224]      [ -3.07207606  -1.35201728   0.35500973]      [  2.63202993  -1.35201728   5.32397834]      [  5.19884868  -1.35201728   1.63816326]      [  7.65721835  -1.35201728   1.13843392]      [  2.48172754  -0.35201728   6.65584512]      [  6.0905911   -0.35201728   1.15552652]      [  8.62497546  -0.35201728   0.30407144]      [ -4.7300089    0.64798272   0.31481496]      [ -3.03274093   0.64798272   0.95337568]      [  2.19653614   0.64798272  10.3675204 ]      [  6.20384058   0.64798272   1.42106077]      [ -4.08636605   1.64798272   0.28875288]      [  2.03344989   1.64798272  13.04648211]      [ -4.11717795   2.64798272   0.39713141]      [  1.93304283   2.64798272  10.41313242]      [ -4.37994815   3.64798272   0.84588643]      [  1.66081408   3.64798272  14.96380955]      [ -4.19024027   4.64798272   0.73216113]      [  1.60252433   4.64798272  14.72419286]      [  6.77837359   4.64798272   0.6186005 ]      [ -4.14362668   5.64798272   0.93673165]      [  1.55372968   5.64798272  12.9421123 ]      [ -4.62223541   6.64798272   0.6510101 ]      [  1.527865     6.64798272  10.80209351]      [  6.86820685   6.64798272   0.82550801]      [ -4.68259732   7.64798272   0.45321369]      [  1.36167494   7.64798272   6.45338514]      [ -5.19205787   8.64798272   0.23935013]      [  1.21003466   8.64798272  10.13528877]      [  7.6689546    8.64798272   0.32421776]      [ -5.36436818   9.64798272   0.79809416]      [  1.26248534   9.64798272   7.67036253]      [  7.35472418   9.64798272   0.92555691]      [ -5.61723652  10.64798272   0.4741007 ]      [  1.23101086  10.64798272   7.97064105]      [ -7.83024735  11.64798272   0.47557318]      [  1.20348982  11.64798272   8.20694816]      [  1.14422758  12.64798272   9.26244889]      [  9.18164464  12.64798272   0.72428381]      [  1.0827069   13.64798272  10.08599118]      [  6.80116007  13.64798272   0.4571425 ]      [  9.384236    13.64798272   0.42399893]      [  1.04053491  14.64798272  10.48370805]      [  9.16197972  14.64798272   0.39930227]      [ -9.85958581  15.64798272   0.39524976]      [  0.9942501   15.64798272   8.39992264]      [  8.07642416  15.64798272   0.61480371]      [  9.55088151  15.64798272   0.54076473]      [ -7.13657331  16.64798272   0.32929172]      [  0.92606211  16.64798272   7.83597033]      [  8.74291069  16.64798272   0.74246827]      [ -7.20022443  17.64798272   0.52555351]      [  0.81344517  17.64798272   6.81654834]      [  8.52844624  17.64798272   0.70543711]      [ -6.97465178  18.64798272   1.04527813]      [  0.61959631  18.64798272  10.33529022]      [  5.733054    18.64798272   1.2309691 ]      [  8.14818453  18.64798272   1.37532423]      [ -6.82823664  19.64798272   2.0314052 ]      [  0.56391636  19.64798272  13.61447357]      [  5.79971126  19.64798272   0.30148347]      [  8.01499476  19.64798272   1.72465327]      [ -6.78504689  20.64798272   2.88657804]      [ -4.79580634  20.64798272   0.36201975]      [  0.548376    20.64798272   7.8414544 ]      [  7.62258506  20.64798272   1.52817905]      [-10.50328534  21.64798272   0.90358671]      [ -6.59976138  21.64798272   2.62980169]      [ -3.71180255  21.64798272   1.27094175]      [  0.5060743   21.64798272  11.06117677]      [  4.51983105  21.64798272   1.74626435]      [  7.50948795  21.64798272   3.46497629]      [ 11.10199877  21.64798272   1.78047269]      [-10.15444935  22.64798272   1.47486166]      [ -6.26274479  22.64798272   4.73707852]      [ -3.45440904  22.64798272   1.72516012]      [  0.52759064  22.64798272  12.58470433]      [  4.22258017  22.64798272   2.63827535]      [  7.03480033  22.64798272   3.506412  ]      [ 10.63560314  22.64798272   3.56076386]      [ -5.95693623  23.64798272   2.97403863]      [ -3.66261423  23.64798272   2.31667236]      [  0.52051366  23.64798272  12.5526344 ]      [  4.21083787  23.64798272   1.95794387]      [  6.82438636  23.64798272   4.77995659]      [ 10.18138299  23.64798272   5.21836205]      [ -9.94629932  24.64798272   0.4074823 ]      [ -5.74101948  24.64798272   2.60992238]      [  0.52987226  24.64798272  10.68846987]      [  6.29981921  24.64798272   3.56204471]      [  9.96431168  24.64798272   2.85079129]      [ -9.64229717  25.64798272   0.4503241 ]      [ -5.579063    25.64798272   0.64475469]      [  0.52053534  25.64798272  10.05046667]      [  5.79167815  25.64798272   0.92797027]      [ 10.05116919  25.64798272   2.52194933]      [ -8.55286247  26.64798272   0.94447148]      [  0.45065604  26.64798272  10.97432823]      [  5.50068393  26.64798272   2.39645232]      [ 10.08992273  26.64798272   2.77716257]      [-16.62381217  27.64798272   0.2021621 ]      [ -9.62146213  27.64798272   0.62245778]      [ -7.66905507  27.64798272   2.84466396]      [  0.38656111  27.64798272  10.74369366]      [  5.76925402  27.64798272   1.13362978]      [  9.83525197  27.64798272   1.18241147]      [-15.64874512  28.64798272   0.18279302]      [ -7.52932494  28.64798272   2.94012191]      [  0.32171219  28.64798272  10.73770466]      [  9.4062684   28.64798272   1.41714298]      [-12.71287717  29.64798272   0.70268073]      [ -7.59473877  29.64798272   2.16183026]      [  0.20748772  29.64798272  12.97312987]      [  3.92952496  29.64798272   1.54987681]      [  9.05148017  29.64798272   2.40563748]      [ 14.96021523  29.64798272   0.55258241]      [-12.14428813  30.64798272   0.36365363]      [ -7.12360666  30.64798272   2.54312163]      [  0.40594038  30.64798272  12.64839117]      [  4.59465757  30.64798272   1.23496581]      [  8.54333134  30.64798272   2.18912857]      [-10.6296531   31.64798272   1.4839259 ]      [ -7.09532763  31.64798272   2.0113838 ]      [  0.37037733  31.64798272  12.2071139 ]      [  3.01253349  31.64798272   3.01591777]      [  4.64523695  31.64798272   3.50267541]      [  8.39369696  31.64798272   2.53195817]      [ -7.07947026  32.64798272   1.01324147]      [  0.39269437  32.64798272   9.67368625]      [  8.58669997  32.64798272   1.00475646]      [ 12.02329114  32.64798272   0.50782399]      [-10.13060786  33.64798272   0.31475653]      [ -7.30360407  33.64798272   0.35065243]      [  0.49556923  33.64798272   9.66608818]      [ -5.37822311  34.64798272   0.38727401]      [  0.4958055   34.64798272   7.5415026 ]      [  6.07719006  34.64798272   0.63012453]      [ -4.64579055  35.64798272   0.39990249]      [  0.46323666  35.64798272   4.60449213]      [  4.72819312  35.64798272   0.98050594]      [ -4.62418372  36.64798272   0.64160709]      [  0.48866236  36.64798272   4.29331656]      [  5.06493722  36.64798272   0.59888608]      [  0.49730481  37.64798272   1.32828464]      [ -1.31849217  38.64798272   0.70780886]      [  1.70966455  38.64798272   0.88052135]      [  0.06305774  39.64798272   0.47366487]      [  2.13639356  39.64798272   0.67971461]      [ -0.84726354  40.64798272   0.63787522]      [  0.55723562  40.64798272   0.62855097]      [  2.22359779  40.64798272   0.33884894]      [  0.77309816  41.64798272   0.4605534 ]      [  0.56144565  42.64798272   0.43678788]]      When plotted, my pointcloud typically looks like this, with color representing the "fitness" of a point (how much its neighbourhood is symmetric according to a certain criterion): ![enter image description here](http://i.stack.imgur.com/4jxIa.png) EDIT: the colored points are obtained as follows:   1. For each horizontal level (points with same Y coordinate) an almost horizontal slice is taken from the back surface model. For each point in this slice, a tangent is drawn and the slice is divided in right side and left side. Then, the difference between the distances of each bilateral pair of points to that line is measured and summed. Then, for each point we have a weight that tells how much the right side and the left side of the slice at that point are equal. The images below illustrate it: Graph showing a slice of the back surface at the lumbar level, where each point receives a value for its asymmetry. The local minima are selected (red points) with a weight (green triangles) proportional to how much it is "better" than its vicinity (how deep its asymmetry value is, locally). ![enter image description here](http://i.stack.imgur.com/7d3Lj.png) I am wondering which fitting/regression method should I use, in order to get:   * Local differentiability (could be approximated);   * Ability to interpolate to arbitrary resolution;   * No overshoot or unpredictable oscillations;   * If possible, outlier removal. I suspect something might be done with least-squares, radial basis functions, RANSAC, etc. but I also suspect I don't have enough knowledge to select a proper method, so that's why I'm asking the experts. EDIT: performing a weighted spline interpolation, I get results similar to the one below, but it is too sensitive to "outliers", and subject to oscillations obviously far away from what should be the "right" central line: ![enter image description here](http://i.stack.imgur.com/xPhkv.png) Any help is much welcome, thanks for reading"


title: Can someone explain what a plim is in regression?

text: I know it has to do with consistency in the coefficients of the explanetory variable but I just don't understand what it is or how we use it. My professor doesn't explain it well, and I cant seem to find a good source that can explain it well. Can someone help?


title: Change in Proportion - Multiple Periods

text: So I have data consisting of a binomial variable observed in multiple time periods and I want to use a statistical test to determine the existence of a change over time. So here's how it looks:               period   number of observations     number of successes    number of fails     1        10                         7                      3     2        30                         22                     8     3        25                         20                     5     4        20                         13                     7     5        4                          2                      2     6        1                          0                      1     7        10                         6                      4      Let's say I want to compare period 5-6-7 with 1-2-3-4. If we group 5-6-7 and 1-2-3-4 together and simply use a $p_{1234}$ and $p_{567}$, where $p_{1234}$ is the success rate in periods 1-2-3-4 and $p_{567}$ is that of periods 5-6-7, we can use the usual $ \frac{p_{1234}-p_{567}}{\sqrt{p_{1234567}(1-p_{1234567}(1/n_{1234}+1/n_{567})}}$ formula. However this does not show that $p_5 < p_{1234}, p_6 < p_{1234}, p_7< p_{1234}$ simultaneously. I sense that I must contruct a suitable F-test but I cannot come up with the exact formulation. How can this be done?


title: What would be an ideal measure for the synchrony of laughter?

text: I am exploring the relationship between the quality of connection two people share and the tendency for their laughter to synchronize. In order to do this, I need some way to quantify the degree of synchrony. The measure that has seemed most promising to me so far is "r" (order) from Kuramoto's model. My understanding is that in using this measure, I am making the assumption that the laughter oscillations are described by a sinusoidal function. My actual data is a series of points, each representing the time at which the amplitude of a laugh-call peaks, so I am wondering if this is an ideal representation.


title: PCA for dimensionality reduction of audio data

text: I am trying to get spectrogram from audio file, using voicebox tool. My input audio file is **audiofile.wav**. After spectrogram, I want to apply PCA to my spectrogram matrix, **b** (to reduce dimensionality). Is it the correct way to apply PCA?               [s,fs]=wavread(audiofile.wav);     [t,f,b]=spgrambw(s,fs,'gmJcw');     coeff=pca(b);     my_features  = b*coeff(:,1:50);      Any suggestions/hint are welcomed.


title: Estimation from two observations

text: Suppose there are two vector signals $x$, $z$. The observer 1 receives a linear version of $x$ plus Gaussian noise. Observer 2 receives a linear sum of both $x$ and $z$ plus Gaussian noise as shown below; $$y_1(t)=A_1(t)x(t)+n_1(t),$$ $$y_2(t)=A_2(t)x(t)+B_2(t)z(t)+n_2(t).$$ But observer 2 is only interested in estimating $z$ so ideally he would like to cancel the whole term $A_2x$. Furthermore $A_2$ is changing slowly over time. Slowly in the sense the changes are correlated over time. Suppose observer 1 can provide his estimates of $x$ say $x_1$ to observer 2. Is there any way for observer 2 to use that information to better estimate $z$? The matrices $A_1,A_2,B_2$ are independent. Edit: Suppose $A_2(t+1)= \sum_{n=0}^{n=t}\alpha_nA_2(n) + z_t(t+1)$ where $z_(t)$ are independent Gaussian noise.


title: Calculate the confidence interval for the mean of a beta distribution

text: Consider a beta distribution for a given set of ratings in [0,1]. After having calculated the mean: $$ \mu = \frac{\alpha}{\alpha+\beta} $$ Is there a way to provide a confidence interval around this mean?


title: What is grid decomposition?

text: I'm reading a paper about DBSCAN implementation on Map-Reduce. In preprocessing step, it wants to divide dataset and it says: > We consider two strategies in this paper for different scenarios. The ﬁrst > one is given the number of computing slots N and its grid decomposition (N = > a ∗ b), it could ﬁnd an approximate division in raw spatial data. I don't get it exactly. Does he want to get coefficients of number of processing unit and use them to divide different dimensions of data? **P.S:** Please help me choose appropriate tags.