title: How to visualise observations of a stochastic process with a parameter change?

text: I've obtained 2840 observations of a stochastic process with large Gaussian noise. The process parameters change at one point around the 1500th observation. I've ensemble-averaged these into sets of 10, 40 and 284 observations and would like to present those three sets to the readers of my report, but I can't figure out a good way of visualising this data. Ideally, the visualisation would let the viewer see:   * the wave morphology   * the impact of noise, i.e. how much one averaging differs from another in the set   * the point where the wave morphology changes I could plot them all in a single graph. This gives a good idea about the impact of noise, but loses the chronological order of the waves. Also, it would only work for the set of 10 before becoming unmanageable. ![enter image description here](http://i.stack.imgur.com/qiL2w.png) I could also plot them as a 3D surface with the X axis showing time in a sample and the Z axis the ordinal number of the observation. One problem is that some waves are seen better than others, making it hard to compare the noise levels and wave morphologies quantitatively from the graph. ![enter image description here](http://i.stack.imgur.com/LFL83.png) How could I best visualise these averaged sets of observations?


title: Measuring temporal aggregation

text: I have a series of 2000 individuals, for each of which I measure the daily duration of a repeating event over the course of a year. That is, for each individual I have a vector like:               id.1020<-c(0,0,0,0,1000,2,15,20.....)       of 364 total length, each element of id.2010[i] representing the duration of the observed event at day[i] I would like some measure of the temporal aggregation of the events duration. That is, individuals that have a similar annual duration (=sum of all daily durations, i.e `sum(id.1020)`, can have their total event duration spread over many days and some over less days. Is there any measure of this "temporal" aggregation that I could use to quantify this "temporal" spread? Thank you very much


title: Is there always a maximizer for any MLE problem?

text: I wonder if there is always a maximizer for any maximum (log-)likelihood estimation problem? In other words, is there some distribution and some of its parameters, for which the MLE problem does not have a maximizer? My question comes from a claim of an engineer that the cost function (likelihood or log-likelihood, I am not sure which was intended) in MLE is always concave and therefore it always has a maximizer. Thanks and regards!


title: Is there a desription in the literature of a Normal hierarchical model with hyperparameters for both the mean and the standard deviation?

text: I'm looking for a comprehensive description of and justification for a Normal hierarchical model where _both_ the means of the groups and the standard deviation are modelled. It is common to find something like the following model in many textbooks (e.g. Gelman et al., p. 288): $$y_{ij} \sim \text{Normal}(\mu_i,\sigma) \\\ \mu_i \sim \text{Normal}(M, S)$$ where $y_{ij}$ is the $i$th datapoint from group $j$ and where non-informative priors are proposed for $M,S$ and $\sigma$. What I'm looking for is an extension of this model where also the standard deviations of each group are modelled and given hyperparameters (and not only a single $\sigma$ is assumed for all groups). That is something like: $$y_{ij} \sim \text{Normal}(\mu_i,\sigma_i) \\\ \mu_i \sim \text{Normal}(M, S) \\\ \sigma_i \sim \text{SomeDistribution}(P_1,P_2,\dots)\\\ $$ but where proposals are given for   * The distribution of the $\sigma_i$s ($\text{SomeDistribution}$ in the model above).   * Non-informative prior distributions for the hyperparameters $M, S$ and the parameters of $\text{SomeDistribution}$. I have not been able to find this in the literature, and my question is: **Where can I find this model described in the literature?** Or alternatively: **What should such a model look like?** ## References Gelman, A., Carlin, J. B., Stern, H. S., Dunson, D. B., Vehtari, A., & Rubin, D. B. (2013). _Bayesian data analysis_. CRC press.


title: 5 Likert Scale and Normal Distribution prior to parametric testing

text: I have data from a 5 point Likert scale item which approximately 1/4 of the responses only account for 1 through 4 of the 1 through 5 available. Do I have to standardize the min/max response categories that is reported in SPSS descriptives for output or will creating z-scores be sufficient? The variables have differing degrees of distributions due to this and the transform functions (e.g., SQRT, LG10) have not solved my problem. My sample (N=162) is the first part of a test/retest for instrument reliability. I have run all the processes thinking that what looked normal and based on criterion (2 or 3 times the SD is sufficient to estimate skewness to identify 'normal distribution' regardless of visual). I'm at the point of frustration and exhaustion...any expert advice?


title: Basic question regarding probability

text: I quote this question from "All of Statistics": > Suppose a gene can be type A or type a. There are three types of people > (called genotypes): AA, Aa, and aa. Let (p, q, r) denote the fraction of > people of each genotype. We assume that everyone contributes one of their > two copies of the gene at random to their children. We also assume that > mates are selected at random. The latter is not realistic however, it is > often reasonable to assume that you do not choose your mate based on whether > they are AA, Aa, or aa. (This would be false if the gene was for eye color > and if people chose mates based on eye color.) Imagine if we pooled > everyone’s genes together. The proportion of A genes is P = p + (q/2) and > the proportion of a genes is Q = r \+ (q/2). A child is AA with probability > P^2, aA with probability 2PQ, and aa with probability Q^2. Thus, the > fraction of A genes in this generation is P^2 + PQ This is probably a basic question. I'm asking: why is A = P^2 + PQ and not P^2 + 2PQ?


title: ARMA model coefficient standard errors

text: I'm writing Python code to use the Kalman State-space approach to estimate ARMA model coefficients using MLE however, I'm not too clear on how to derive the coefficient estimates standard errors from the process. I know you have to use the information matrix, invert it to get the variance- covariance matrix and you can derive it from there, but how to derive the information matrix is the challenge I'm currently faced with (from a coding standpoint). Any links to free papers, tutorials, equations, or references to code is what would be very helpful to me. Thanks!


title: Organizing a classification tree (in rpart) into a set of rules?

text: Is there a way that once a complex classification tree is constructed using rpart (in R), to organize the decision rules produced for each class? So instead of getting one huge tree, we get a set of rules for each of the classes? (if so, how?) Here is a simple code example to show examples on:               fit <- rpart(Kyphosis ~ Age + Number + Start, data=kyphosis)      Thanks.


title: Multivariate mutual independence test

text: I have a series of observations $\\{X_i,Y_i,Z_i\\}$ for random variable $X$, $Y$ and $Z$. Now I want to test if $X$, $Y$ and $Z$ are mutual independent, can anyone help me? Some pairwise independence tests like this one http://www.gatsby.ucl.ac.uk/~gretton/indepTestFiles/indep.htm are available, but I haven't found such tests for multivariate mutual independence yet? Thanks in advance.


title: Panel data model for exports and exchange rates

text: Suppose I have 4 years worth of _monthly_ panel data on:   1. exports of widgets $y$ from home country to 12 different nations (in US dollars)   2. nominal exchange rates $x$ for those 12 countries (in US dollars per unit of local currency)   3. domestic consumption of widgets $z$ (in US dollars)   4. nominal exchange rate $r$ of home country (in US dollars per unit of home currency) I am primarily interested in the effect of exchange rate volatility on exports. I estimated the following seasonally-differenced equation using OLS with standard errors clustered on country and time: $$ \ln y_{i,t} - \ln y_{i,t-52} = \alpha + \beta \cdot(\ln x_{i,t}-\ln x_{i,t-52})+\gamma \cdot(\ln z_{i,t}-\ln z_{i,t-52}) +\eta \cdot(\ln r_{i,t}-\ln r_{i,t-52}) + (\varepsilon_{i,t}-\varepsilon_{i,t-52})$$ The subscript $i$ indexes countries and $t$ indexes months. Exports are all strictly positive, so taking logs is not a problem. Two Questions:   * My estimates are somewhat sensitive to whether I control for own exchange rate $r$, and $\hat \beta \approx -\hat \eta$. Does $r$ belong in the equation?   * Are the any robustness checks that I should consider? Are there any problems that jump out? I am not familiar with this literature, so any references would be most welcome.


title: Multilevel - focus on higher level

text: I have a data structure similar to firms (the higher level) with individuals within them (the lower level). Each firm can have any number of individuals present. I have a reasonable set of predictors for both levels but my aim is to produce estimations for a dependent variable at the firm level. I have no reason to believe that the firms are not independent from each other. It would seem logical that the individuals could be correlated with respect to firm: certain individuals are more likely to cluster together at certain firms. My question is - am I right in thinking that since my focus here is on the higher level, I can simply aggregate predictors held about individuals to the firm-level and apply them as predictors against the firm? My understanding is that the dependence between individuals can be captured within the aggregation but since firms are independent from one another, this will not bias the standard errors.


title: how to keep covariates when creating a sequence object from a SPELL data

text: I have a SPELL data in which there are many covariates together with the basic variables (i.e., "id", "begin", "end", and "status"). When I define the SPELL data as a sequence object by "seqdef", it seems that those covariates disappear in the sequence object. Due to this, I cannot make graphs to be grouped by their levels. Do you have any suggestions?


title: Triple nested cross validation

text: I have read several very informative posts including the link about the nested/double cross validation, which can determine (sub)optimal hyperparameter values as well as make an unbiased estimate of the generalized performance. I wonder if some sort of "triple" cross validation would make sense for accomplishing the following:   1. selecting hyperparameter values,   2. selecting best features/predictors/attributes (feature selection),   3. getting an unbiased estimate of the generalized performance. The workaround would probably be to treat features as parameters, make a bigger grid (with combinations of both the features and parameters) and apply nested/double cross validation. Are there any "standard" approaches for accomplishing all three aforementioned tasks with cross validation?


title: DIDM Matching with Difference-in-difference using cross-sectional data

text: I need to know how can I combine in stata the DID (difference-in-difference) model with the PSM (Propensity score matching) with cross-sectional data. I had searched and there are tones of answers related to each of those models, but together I haven't been able to find them. Could someone please help me!! ??


title: Non parameteric tests for 2*2 mixed design

text: I have 2*2 mixed design, where the two factors are _time_ and _group_. The _time_ levels are _pre-_ and _post-intervention_ , and the _group_ levels are _control_ and _intervention_. If I were using a parameteric test, I think a mixed ANOVA would be appropriate, but as the data is ordinal, I thought a non-parametric test would be more suitable. Am I right in thinking that there isn't a non-parametric equivalent of a mixed ANOVA? Is using a Wilcoxon Signed Ranks test (for a within-subjects comparison of the _pre_ and _post_ ) in conjunction with a Mann-Whitney U test (for a between-subjects comparison of the _control_ and _intervention_ ) an appropriate substitute? If so, would I use Bonferonni's correction when looking at the p-levels (because I'm running two tests)? Many thanks.


title: Assumptions of the Ordered Probit model

text: What are the assumptions of an ordered probit model that must be met? What are the tests to check these? Many thanks in advance


title: Flexible and inflexible models in machine learning

text: I came across a simple question on comparing flexible models (i.e. splines) vs. inflexible models (e.g. linear regression) under different scenarios. The question is: In general, do we expect the performance of a flexible statistical learning method to perform better or worse than an inflexible method when:   1. The number of predictors $p$ is extremely large, and the number of observations $n$ is small?    2. The variance of the error terms, i.e. $σ^2 = \text{Var}(e)$, is extremely high? I think for (1), when $n$ is small, inflexible models are better (not sure). For (2), I don't know which model is (relatively) better.


title: Approaching the limit in the .632 rule when n is unknown

text: This question relates to the ".632 rule" asked about in an earlier question here. My question is, assuming we have a process that approximates the random- sample-with-replacement for large (unknown) n contemplated by the rule, how can we measure our progress along the asymptotic curve toward $1/e?$ I had an inchoate idea, basically that we keep track of the frequency of previously noted marbles (so to speak) and use either the change in frequency over a given number of samples or use the absolute number _ab initio_ as a measure of progress. The bottom line is a pragmatic determination of a point of diminishing return. Roughly speaking, if I am a biologist sampling a pond's population of rotifers, how do I know when I have seen about all I am likely to see? Thanks.


title: Inverse matrix for contrast coding

text: I'm trying to understand how "user defined contrast coding" works. My question refers to the example from http://www.ats.ucla.edu/stat/r/library/contrast_coding.htm#User:               #initial contrast matrix including the constant term     mat = matrix(c(1/4, 1/4, 1/4, 1/4, 1, 0, -1, 0, -1/2, 1, 0, -1/2, -1/2, -1/2, 1/2, 1/2), ncol = 4)     mat          [,1] [,2] [,3] [,4]     [1,] 0.25    1 -0.5 -0.5     [2,] 0.25    0  1.0 -0.5     [3,] 0.25   -1  0.0  0.5     [4,] 0.25    0 -0.5  0.5          mymat = solve(t(mat))     mymat          [,1] [,2] [,3] [,4]     [1,]    1 -0.5   -1 -1.5     [2,]    1  0.5    1  0.5     [3,]    1 -1.5   -1 -1.5     [4,]    1  1.5    1  2.5          #remove the intercept (constant) term     my.contrasts<-mymat[,2:4]     contrasts(hsb2$race.f) = my.contrasts      **Question:**   * Why is it necessary to calculate the inverse of the transposed matrix? All other examples on the page doesn't use matrix algebra (e.g "Dummy Coding").


title: Renyi divergence identity

text: I'm reading the paper, T. van Erven and P. Harremoës, Rényi Divergence and Kullback-Leibler Divergence, _arXiv 1206.2459_ on the Renyi divergence, and I'm trying to make sense of "Example 1". I think that I'm just getting messed up on the terminology, so I'm not able to derive that the identity is true. Am I missing some assuptions/definitions/notational conventions so that this doesn't make sense? Is there a different way of describing/deriving their "Example 1"? > **Example 1** : Let $Q$ be a probability distribution and $A$ a set with > positive probability.   >  Let $P$ be the conditional distribution of $Q$ given $A$.   >  Then $$ D_\alpha( P \Vert Q ) = - \ln Q(A). $$ Some features of their notation (from earlier in the paper): $$ P=(p_1, p_2, \dots, p_n) \\\ Q=(q_1, q_2, \dots, q_n) $$ refer to (discrete) probability distributions. The Renyi divergence itself is expressed as: $$ D_\alpha ( P \Vert Q ) = \frac{1}{1-\alpha} \ln \sum_{i=1}^{n} p_i^{\alpha} q_i^{1-\alpha} $$ I don't understand the notion of "the conditional distribution of $Q$ given $A$", to me, this reads as "the probability distribution of a probability distribution [function]". I'm also unclear on the interpretation the result: the result should be the number, but what's written (looks like it) is the log of a set (note the definition of $A$).


title: What is quantile normalization?

text: Is anyone familiar with quantile normalization? As I understand it's a way to "polish" a set of datasets so that they have the same "statistical properties" (according to Wikipedia)? I wonder if anyone knows what are the statistical properties that they refer to when they say > quantile normalization is a technique for making two distributions identical > in statistical properties. Are those statistical properties the mean, median, or standard deviation?


title: Predictive performance depends more on expertise of data analyst than on method?

text: I've come across a rumour that some study showed that the performance of predictive models depends more on the expertise of the data analyst with the chosen method than on the choice of the method.   In other words, the claim is that it is more important that the data analyst is familiar with the chosen method than how "appropriate" the method would seem for the problem from a more theoretical standpoint. This was mentioned in the context of chemometrics, which involves typically problems of many variates (100s - 1000s), multiple collinearity, and of course, too few samples. Prediction may have been classification or regression. My personal experience suggests that this is _plausible_ , but a study was mentioned (I asked the person who mentioned that by email after a quick but unsuccessful search, but never received any answer). However, also with a more elaborate search, I was not able to track down any papers. Is anyone aware of such findings? If not, what does the personal experience of Big Guys here say?


title: Finding kind of dependency

text: I have an interesting (to me) real-life problem but only a very basic understanding of statistics. I would appreciate a pointer to literature as I understand actually solving the problem is my responsibility :) The model: There are M * N coins categorized by two nominal attributes A and B with members' counts M and N correspondingly - so that for every combination of attributes' values there is exactly one coin - basically, an M * N matrix of coins. It is hypothesized that bias of every coin is a function of A and B - let's denote it f(A, B). The nature if this function is not exactly known. The data: Every coin is flipped K times, where K for different coins are not necessarily equal (some may be even zero), so K is a matrix M * N. The number of heads and tails is H and T (matrices as well). Ultimately, I need to estimate biases of all coins, including those with zero data (which I guess is possible due to the fact that biases are mutually dependent), but before getting there (I think) I need to establish a good enough hypothesis for f(A, B). What would be a recommended methodology to find the f(A, B) given matrices of H and T (from multiple experiments, if needed)? I would imagine some kind of exploratory analysis would be needed to form a few hypotheses about possible kinds of functions, and then Bayesian inference could be used to find the specific parameters of the function? If it helps, I am willing to assume f(A, B) can be factored as g(c(A), d(B)), where c and d are scalar functions. Thanks in advance! PS: I do not control K directly, they are in fact gradually accumulated as a result of a separate stochastic process. I guess it makes sense to solve the current problem before bringing _that_ process into scope, so for know I treat K as fixed.


title: How do I measure how well my data fits a prediction

text: Given a set of data like this:               x     y1            y2     50    102.50        16     100   408.60        64     200   1,871.30      256     400   8,824.40      1024     800   49,695.60     4096      My hypothesis is that when `x` doubles, you should find that `y` quadruples. I believe this is equivalent to saying that `x` is proportional to the square of `y`. * * * We can see that the `y1` data (from my actual problem) are close to this hypothesis, but not in complete agreement; by comparison the `y2` data (which I have contrived) match my hypothesis exactly. What I am looking for is some quantifiable measure of how much my data agrees with my hypothesis, and the correct terminology to discuss this type of test. I would expect that `y2` would be fully agree, and `y1` would provide some distance metric.


title: ARIMA Model building

text: Can someone help me with this please: it's revision and I am extremely stuck with what to do and it's holding me back! A previous exam question asks: >   * explain the following ARIMA printouts and what is their contribution to > selecting the final ARIMA model >   * what would you do next in the stage of model identification and > selection? >   * For model 1, write out the mathematical equation of the built model > (additive). > Given this data:               Trendcycle Model 1:     Normalised BIC = -0.743     Ljung Box Q = Stats (9.808), sig (0.832)     ARIMA model parameters=     Constant = estimate (31.812), t (21.695), sig (0.00)      AR Lag 1 = estimate (0.860), t (11.201), sig (0.00)     Lag 2 = estimate (0.880), t (11.440), sig (0.00)     Lag 3 = estimate (-0.780), t (19.071), sig (0.00)          Trendcycle Model 2:     Normalised BIC = -0.910     Ljung Box Q = Stats (11.544), sig (0.713)     ARIMA model parameters=     Constant = estimate (-0.034), t (-0.174), sig (0.863)      AR Lag 1 = estimate (-0.154), t (-1.108), sig (0.273)     Lag 2 = estimate (0.778), t (8.418), sig (0.00)     Lag 3 = estimate (0.048), t (0.326), sig (0.746)      Thanks everyone. Help needed :D


title: Modeling measurement error (reliability?) using repeated measurements

text: I have a measurement protocol that seeks to characterise the visual sensitivity of a subject. I would like to assess how reliable the measurements from this procedure are. I have measured a small number of subjects, with two repeated measurements per day for a number of days (not all subjects did the same number of days). The data are similar in form to this (code for R):               # generate some example data:     dat <- expand.grid(subject=1:4,rep=c(1,2),day=seq(0,28,by=4))          # extra subject with more days:     d <- expand.grid(subject=5,rep=c(1,2),day=seq(0,36,by=4))     dat <- rbind(dat,d)          dat$subject <- factor(dat$subject)          dat$y <- rnorm(nrow(dat),mean=3,sd=0.25)      which produces some data whose summary looks like:                subject      rep           day              y              1:16    Min.   :1.0   Min.   : 0.00   Min.   :2.455        2:16    1st Qu.:1.0   1st Qu.: 8.00   1st Qu.:2.817        3:16    Median :1.5   Median :16.00   Median :2.970        4:16    Mean   :1.5   Mean   :14.95   Mean   :2.978        5:20    3rd Qu.:2.0   3rd Qu.:24.00   3rd Qu.:3.141                Max.   :2.0   Max.   :36.00   Max.   :3.617       The measurement value is ratio-scaled. Day is relative to the other days tested within a subject (i.e. not all subjects' "day N" corresponds to the same actual day). One simple measure of consistency is to correlate the scores within a day - a simple test-retest measurement. However, this doesn't usefully use the information that the same subject repeated the measurement across days. Another measure that I considered was an autocorrelation of the measurements as a function of time, but this is extremely noisy because not all subjects did the same number of days, so it must be done for each subject. I also modelled the data using a linear model with "days" and "rep" as continuous predictors, in a Bayesian hierarchical framework. This shows that the slope for day and the slope for repetition within a day are not credibly different from zero. This indicates that the test scores are fairly stable across days and repetitions. However, one thing I would like to be able to do is to directly compare the variance between subjects with the variance within a day or across days within subjects. That is, to test whether the variability of measurement error + within subject fluctuations is smaller than the variability in normal visual sensitivity across individuals. Does anyone have any suggestions for how to model this?


title: How to interpret the output of a multinomial classification model in R package gbm

text: After running a gradient boosted model with `n` data points using multinomial regression where the response variable (a factor, as required by the gbm function) has `k` levels with R package gbm, I see that the predictions are output as as a vector of length `n*k`. Predicted responses are from:               probs.var.multinom <- predict.gbm(gbm.model.multinom, test.data, best.iter.gbm, type="response")      Note that this is different from the output of a logistic (distribution = "bernoulli") model, where the results are a vector the same length as the number of cases. How should this be interpreted? Specifically, how can I link the response vector back to the input data set to evaluate the classification?


title: On cross-validation schemes for "rectangular" samples

text: Consider this example. Suppose that for any pair $(x, y)$ of bacterial strain $x$ and (candidate) anti-bacterial agent $y$, we can experimentally determine some measure $f(x, y)$ (say, the concentration of $y$ needed to kill off 50% of a culture of $x$ in some standardized assay). Now we want to devise a statistical model $M$ to predict $f(x, y)$ on the basis of some features of $x$ and $y$. We decide on sets of features $V$ and $W$ for bacterial strains and anti-bacterials, respectively. Therefore, each pair $(x, y)$ of strains and agents gets assigned a composite feature vector $(\mathbf{v}(x), \mathbf{w}(y))$. The goal is to optimize the estimates $M(\mathbf{v}(x), \mathbf{w}(y))$ of $f(x, y)$. To train the model, our "gold standard" dataset consists of the experimentally measured $f(x, y)$ values for all $(x, y)$ pairs in the Cartesian product $X \times Y$ of some panel $X$ of strains and panel $Y$ of agents. (This, BTW, is the "rectangular" sample alluded to in the post's title.) To prevent over-fitting, we want to use some form of cross-validation. My question concerns the protocol for selecting the subsets $T$ and $V$ of $X \times Y$ to use for training and validation, respectively. My intuition is that $T$ should be in the form of a Cartesian product $X_t \times Y_t$, where $X_t$ and $Y_t$ are randomly chosen _proper_ subsets of $X$ and $Y$, respectively; and that $V$ should be of the form $X_v \times Y_v$, where $X_v = X - X_t$ and $Y_v = Y - Y_t$. Of course, this scheme does not use all of the data. An alternative scheme would be to make the training set $T$ as a random subset of $X \times Y$, and take $V = (X \times Y) - T$. This second scheme, however, seems wrong to me, because such training set $T$ would contain too much information about all of the strains in $X$ and the agents in $Y$ to make $V$ a fair challenge to the model. Is this the right intuition? If so, what's the rigorous way to show it?


title: Data analysis on multiuser environment

text: I have following environment that I'd like to analyze. My company is divided in two teams:   * Insert Team (composed by 150 persons)   * Quality Check Team (composed by 10 persons) Suppose that all people of first team compile a form with 10 text boxes, for example: name, surname, address, category, studies, preferences, hobbies, etc., and the second team (quality) analyzes a variable percentage of the data which are inserted. Now I'd like to know which method could I apply to find good or poor members of first team, as analyzed by second team, **captured using historical feedback of all of the analysis made before**. An example to understand better: If the first team inserts 10 data points for the users, then we'll have 1500 data points to analyze. Suppose we define 10% of those data to pass to the second team. After their analysis and feedback, I'd like to use it to influence the next 10% of data capturing, for example, to focus attention on users that are less precise or on text box where there are more errors (which is the approach of major companies like Facebook or Twitter).


title: Why is power analysis with logistic regression so liberal compared to chi squared?

text: I'm conducting a power analysis to derive the required sample size for a study - basically compared exposed / non-exposed with 30-day mortality as outcome. I'll check for crude mortality rates with chi-square, but also use logistic regression with probable confounders. When I run a power analysis - power 0.8, significance level 0.05, effect size 0.15 and estimated 10 confounders I get that I'd need only n=117 which seem quite small. comparing with chi-square - it suggest that I'd need 350. I'm using R and `pwr`:               pwr.f2.test(u=10, v=NULL, f2=0.15, sig.level=0.05, power=0.8)     pwr.chisq.test(w=0.15, N=NULL, df=1 , sig.level=0.05, power=0.8 )      Is this predictable or am I misusing this?


title: Understanding random forest, gini, and KS

text: I'm a beginner machine learning user, doing my first predictive model using random forest. I have some questions regarding the way to measure how good a model is (Gini area from roc curve, and KS), and some about random forest algorithm. Thanks a lot for your help. 1) Which metric is better to compare models? Gini or KS? For example, I have two versions of my model, one in which gini is greater than the other but lower KS. 2) Once I run my model, I can see how significant are each variable I created using VarImpPlot in R. All of them "explain" something, but I struggle to know "when my model is optimal". For example, some times I remove one of the variables and it improves, or sometimes I remove not one but a combination of two, and it improves... but some times it gets worse. Is there any way to know when I reached an "optimal combination" given my set of variables? Thanks a lot


title: Consistency of unbiased estimator of error term variance in Multiple regression

text: Let $Y=X\beta+\epsilon$. We know that $\frac{e'e}{n-k}$ is an unbiased estimator of $Var(\epsilon)$, where $e$ is the vector of residuals, and $\epsilon$ is multivariate normal distributed in this model. How can we prove, if it is, that this estimator converges to the right value? Or even if it's consistent? Any help would be appreciated.


title: Why compute a Q(s,a) state-action function when the model is unknown in reinforcement learning?

text: Sutton and Barto state: (p.117, "Reinforcement Learning") "With a model (e.g. a full DP or MDP), state values (e.g. V(s) are sufficient to determine a policy- just look ahead to the next step and choose the best combination of reward and state. Without a model, however, state values alone are not sufficient. One must explicitly estimate the value of each action for the values - e.g Q(s,a) to be useful in suggesting a policy." What is an example that makes this obvious? Why can't one just use a form of Bellman's equation in V(s) to determine the best policy via e.g. value iteration. without a state model?


title: Dealing with Euclidean distance and dimension independence

text: I'm not very well informed in terms of distance measures. What sort of distance measure would I use if I know that the various dimensions are not independent of each other? Because I think that euclidean distance assumes that the dimensions are independent.


title: Confidence interval on entropy estimation

text: Over a sample (of about a hundred observations), I have estimated the "empirical" entropy of a sequence of identically independently distributed random variables which have 5 possible outcomes. I was curious to know if there was a possibility to estimate in some way a confidence interval over this "empirical" entropy calculation (i.e. where true probabilities are replaced by their observed frequencies) ? Best regards


title: How to calculate and evaluate estimators in decomposing interaction models to subgroup models?

text: Suppose we have a simple interaction model (ignoring error terms for simplicity) Y = b0 + b1X + b2Z + b3X*Z in which Z is the dummy variable. The model gives values and s.e. of these coefficients. Now we "decompose" this model into sub-group models: For Z=0: Y = b0 + b1X For Z=1: Y = (b0 + b2) + (b1+b3)X We can calculate coefficient values of sub-group models from the interaction model. My question is: Can we calculate s.e. of coefficients of sub-group models from the interaction model in some way? Let's say that we have to do it manually, i.e., without a program. Thanks.


title: Stat test to compare two sets of binomial data

text: I am doing my undergraduate project and have tested a drug on cells for a responce. The data I have collected shows that the cells either produced a responce, or they did not. Binomial data?               Control 0ug/mL      22  3053     Control2 (DMSO)     15  2807     5ug/mL rapamycin    12  3070     10ug/mL rapamycin   7   2806      The first number is the number of cells that produced the responce, and the second number is the total number of cells tested. I'm trying to work out the test I can use to see if the results between the groups are significant. I was think a Fisher's exact test, but I wasn't sure if that was just for 2x2 tables. Would a chi sq test be sufficient?


title: How to simulate a hidden Markov chain?

text: I want to simulate data from a 3-state hidden Markov chain with a known matrix of transition probabilities. Each state corresponds to a bivariate data with known marginals that the dependence between them is modeled by a copula. What is the best algorithm? Is there any function in `R` for such work?


title: Multilevel modeling of experimental data with repeated measures in the same condition of factorial experiment

text: I'm struggling with modeling some experimental data using the `lme4` package in R, and would appreciate input. My experimental design is as follows: subjects entering the experiment answer a screener question, which is used to randomly assign them to a between- subjects condition (variable=`rank`, which has 2-levels (0/1)). They then complete a distractor task, and make two choices (`choice` being the within- subjects dependent variable). At each choice, the stimuli are RANDOMLY ASSIGNED and crossed by two factors (`msg` has 3-levels ("norm", "no norm" and "provincial"), and `cost` has 2-levels (0/1)). **Because participants were randomly assigned to both`cost` and `msg` at two points in time, theoretically they could have the same combination of `cost` and `msg` at both points in time.** My life would be easiest if I could use a wonderful package such as `ezANOVA`, but my data won't allow me to do this because every individual doesn't have EVERY combination of the two within-subjects variables see here. So, I'm in the less-familiar territory of mixed models. My hypothesis argues that there should be a three-way interaction between `rank`, `msg`, and `cost`. Thus, a simple version of my model might be:               m1 <- glmer(choice ~ msg*cost*rank + (1|id), data=df, family="binomial")      But, I've also seen it suggested on this site that my random effects should be modeled as follows in order to evaluate the interaction between the factor and subjects:               m1 = lmer(choice ~ msg*cost*rank + (rank|id) + (cost|id), data=df)      The latter model has me straying into unfamilar territory, so I'd appreciate any advice about how to model this data in a way that is (1) simple, but (2) appropriate. Sample data below:               > dput(df[1:700,2:6])     structure(list(time = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1), choice = c(1, 1,      1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0,      1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1,      1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0,      0, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0,      0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1,      1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1,      1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1,      1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1,      1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1,      1, 1, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1,      1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 1, 1, 1, 0, 1,      1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 0,      1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0,      1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1,      0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1,      1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 1, 1, 0, 1, 1, 0,      0, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 0, 1,      1, 1, 1, 1, 0, 0, 1, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 1, 1, 1,      1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 0, 1, 0,      1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 1,      0, 0, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 1, 1,      1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1,      1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1,      1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0, 1, 0, 1, 0,      0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1,      1, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1,      1, 1, 0, 1, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0,      1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0,      1, 1, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0,      1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1,      1, 0, 0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1,      1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 1, 1, 1, 1, 0, 1,      0, 1, 0, 1, 0), msg = structure(c(3L, 1L, 1L, 2L, 3L, 1L, 3L,      3L, 3L, 2L, 2L, 3L, 2L, 3L, 2L, 3L, 1L, 3L, 2L, 3L, 1L, 1L, 3L,      1L, 2L, 3L, 3L, 2L, 2L, 1L, 1L, 1L, 1L, 3L, 2L, 3L, 3L, 2L, 3L,      3L, 1L, 3L, 1L, 2L, 2L, 3L, 2L, 3L, 3L, 2L, 3L, 3L, 1L, 2L, 3L,      3L, 1L, 2L, 3L, 2L, 3L, 3L, 1L, 1L, 3L, 1L, 3L, 2L, 1L, 3L, 2L,      3L, 3L, 2L, 2L, 1L, 2L, 2L, 2L, 3L, 3L, 2L, 2L, 2L, 2L, 2L, 1L,      3L, 2L, 3L, 2L, 3L, 3L, 3L, 1L, 1L, 3L, 1L, 2L, 2L, 3L, 2L, 3L,      3L, 2L, 2L, 3L, 2L, 1L, 2L, 1L, 3L, 2L, 2L, 1L, 3L, 3L, 2L, 3L,      3L, 3L, 2L, 2L, 1L, 2L, 3L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 2L, 2L,      1L, 1L, 3L, 1L, 3L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L,      2L, 1L, 2L, 1L, 1L, 3L, 3L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 3L,      3L, 1L, 2L, 3L, 1L, 1L, 3L, 2L, 2L, 3L, 3L, 1L, 1L, 1L, 2L, 1L,      2L, 3L, 3L, 2L, 1L, 2L, 3L, 1L, 2L, 2L, 1L, 3L, 3L, 1L, 1L, 1L,      3L, 2L, 3L, 1L, 2L, 2L, 3L, 2L, 1L, 3L, 1L, 2L, 2L, 3L, 3L, 2L,      1L, 3L, 3L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 3L, 2L, 2L, 1L, 1L, 1L,      1L, 1L, 3L, 2L, 2L, 1L, 2L, 2L, 3L, 1L, 1L, 2L, 3L, 2L, 3L, 3L,      3L, 3L, 1L, 3L, 2L, 1L, 2L, 3L, 1L, 1L, 2L, 3L, 3L, 2L, 1L, 1L,      2L, 1L, 2L, 3L, 3L, 3L, 1L, 2L, 2L, 3L, 1L, 3L, 1L, 3L, 3L, 1L,      1L, 3L, 1L, 3L, 1L, 1L, 3L, 1L, 1L, 2L, 2L, 3L, 2L, 3L, 2L, 3L,      1L, 2L, 2L, 1L, 2L, 3L, 3L, 3L, 2L, 3L, 1L, 3L, 1L, 2L, 3L, 1L,      1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 3L, 2L, 3L, 2L, 1L, 3L, 3L,      3L, 1L, 3L, 3L, 1L, 3L, 1L, 3L, 2L, 1L, 3L, 1L, 1L, 2L, 2L, 3L,      1L, 3L, 3L, 3L, 1L, 3L, 1L, 3L, 2L, 1L, 2L, 1L, 1L, 2L, 3L, 3L,      3L, 3L, 3L, 2L, 2L, 3L, 2L, 1L, 1L, 2L, 1L, 2L, 3L, 3L, 1L, 1L,      1L, 3L, 3L, 1L, 1L, 3L, 1L, 2L, 3L, 1L, 1L, 3L, 2L, 1L, 3L, 3L,      3L, 1L, 3L, 1L, 3L, 1L, 3L, 2L, 3L, 2L, 2L, 1L, 2L, 1L, 3L, 2L,      2L, 3L, 1L, 2L, 1L, 3L, 3L, 3L, 3L, 2L, 2L, 3L, 3L, 3L, 1L, 1L,      1L, 1L, 3L, 1L, 3L, 1L, 1L, 1L, 3L, 3L, 2L, 1L, 2L, 3L, 2L, 3L,      2L, 3L, 1L, 3L, 1L, 3L, 1L, 1L, 3L, 1L, 3L, 2L, 3L, 2L, 1L, 1L,      1L, 2L, 2L, 1L, 2L, 1L, 3L, 1L, 2L, 2L, 2L, 1L, 1L, 3L, 3L, 1L,      1L, 3L, 1L, 3L, 3L, 3L, 3L, 2L, 1L, 1L, 3L, 2L, 1L, 3L, 2L, 2L,      1L, 1L, 3L, 3L, 3L, 1L, 1L, 2L, 1L, 1L, 3L, 1L, 2L, 2L, 3L, 1L,      2L, 2L, 3L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 3L, 1L, 3L,      1L, 1L, 3L, 3L, 3L, 2L, 3L, 2L, 2L, 3L, 1L, 2L, 3L, 3L, 2L, 2L,      1L, 2L, 3L, 2L, 1L, 2L, 3L, 3L, 2L, 2L, 2L, 3L, 1L, 2L, 3L, 2L,      3L, 3L, 3L, 3L, 1L, 2L, 3L, 2L, 3L, 3L, 3L, 2L, 3L, 1L, 1L, 2L,      2L, 3L, 3L, 1L, 1L, 3L, 1L, 1L, 1L, 2L, 3L, 2L, 3L, 1L, 2L, 3L,      2L, 2L, 3L, 3L, 3L, 1L, 3L, 2L, 2L, 1L, 2L, 3L, 1L, 3L, 2L, 3L,      1L, 2L, 3L, 2L, 1L, 2L, 3L, 3L, 3L, 3L, 2L, 3L, 2L, 3L, 3L, 2L,      1L, 2L, 3L, 2L, 1L, 2L, 1L, 3L, 1L, 3L, 1L, 1L, 2L, 1L, 3L, 1L,      2L, 2L, 3L, 2L, 2L, 3L, 3L, 1L, 3L, 1L, 3L, 2L, 1L, 3L, 1L, 2L,      1L, 3L, 1L, 3L, 1L, 1L, 2L, 1L, 3L, 2L, 3L, 3L, 3L, 1L, 3L, 3L,      2L, 2L, 1L, 2L, 2L, 2L, 2L, 3L, 1L, 1L, 1L, 2L, 3L, 2L, 3L, 1L,      3L, 3L, 3L, 1L, 3L, 3L, 2L, 2L, 3L, 2L, 2L, 3L, 2L, 1L, 2L, 3L,      2L, 2L, 1L, 3L, 2L), .Label = c("No norm", "Norm", "Provincial"     ), class = "factor"), cost = structure(c(1L, 2L, 1L, 1L,      1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L,      2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L,      2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 2L,      1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L,      2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L,      1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L,      2L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L,      1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L,      2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L,      1L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L,      2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L,      1L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L,      1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L,      1L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L,      2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 1L,      2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L,      1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 2L,      1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L,      1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L,      1L, 1L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,      1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 1L,      2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 2L,      1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L,      1L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L,      1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L,      2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 2L,      1L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 2L,      1L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L,      2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L,      1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L,      2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 1L,      1L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L,      2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L,      2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 1L,      1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L,      1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L,      2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 1L,      1L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L,      2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 1L,      2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 1L,      1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L,      1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 1L,      2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L,      1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L), .Label = c("0", "1"), class = "factor"),          rank = structure(c(1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L,          2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L,          2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 2L,          1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L,          1L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L,          1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L,          1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L,          1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 2L,          2L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L,          2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L,          2L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L,          2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L,          2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L,          2L, 1L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L,          2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L,          2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 1L,          1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L,          1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L,          2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 2L,          2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L,          1L, 2L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L,          2L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L,          1L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L,          1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,          2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L,          1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L,          1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L,          2L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 1L,          2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L,          2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L,          1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 2L,          2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L,          2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L,          2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L,          1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L,          1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L,          1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L,          1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L,          2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L,          1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L,          1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L,          2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L,          2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,          1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L,          1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L,          2L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L,          1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L,          1L), .Label = c("0", "1"), class = "factor")), .Names = c("time",      "choice", "msg", "cost", "rank"), row.names = c("1.1",      "2.1", "3.1", "4.1", "5.1", "6.1", "7.1", "8.1", "9.1", "10.1",      "11.1", "12.1", "13.1", "14.1", "15.1", "16.1", "17.1", "18.1",      "19.1", "20.1", "21.1", "22.1", "23.1", "24.1", "25.1", "26.1",      "27.1", "28.1", "29.1", "30.1", "31.1", "32.1", "33.1", "34.1",      "35.1", "36.1", "37.1", "38.1", "39.1", "40.1", "41.1", "42.1",      "43.1", "44.1", "45.1", "46.1", "47.1", "48.1", "49.1", "50.1",      "51.1", "52.1", "53.1", "54.1", "55.1", "56.1", "57.1", "58.1",      "59.1", "60.1", "61.1", "62.1", "63.1", "64.1", "65.1", "66.1",      "67.1", "68.1", "69.1", "70.1", "71.1", "72.1", "73.1", "74.1",      "75.1", "76.1", "77.1", "78.1", "79.1", "80.1", "81.1", "82.1",      "83.1", "84.1", "85.1", "86.1", "87.1", "88.1", "89.1", "90.1",      "91.1", "92.1", "93.1", "94.1", "95.1", "96.1", "97.1", "98.1",      "99.1", "100.1", "101.1", "102.1", "103.1", "104.1", "105.1",      "106.1", "107.1", "108.1", "109.1", "110.1", "111.1", "112.1",      "113.1", "114.1", "115.1", "116.1", "117.1", "118.1", "119.1",      "120.1", "121.1", "122.1", "123.1", "124.1", "125.1", "126.1",      "127.1", "128.1", "129.1", "130.1", "131.1", "132.1", "133.1",      "134.1", "135.1", "136.1", "137.1", "138.1", "139.1", "140.1",      "141.1", "142.1", "143.1", "144.1", "145.1", "146.1", "147.1",      "148.1", "149.1", "150.1", "151.1", "152.1", "153.1", "154.1",      "155.1", "156.1", "157.1", "158.1", "159.1", "160.1", "161.1",      "162.1", "163.1", "164.1", "165.1", "166.1", "167.1", "168.1",      "169.1", "170.1", "171.1", "172.1", "173.1", "174.1", "175.1",      "176.1", "177.1", "178.1", "179.1", "180.1", "181.1", "182.1",      "183.1", "184.1", "185.1", "186.1", "187.1", "188.1", "189.1",      "190.1", "191.1", "192.1", "193.1", "194.1", "195.1", "196.1",      "197.1", "198.1", "199.1", "200.1", "201.1", "202.1", "203.1",      "204.1", "205.1", "206.1", "207.1", "208.1", "209.1", "210.1",      "211.1", "212.1", "213.1", "214.1", "215.1", "216.1", "217.1",      "218.1", "219.1", "220.1", "221.1", "222.1", "223.1", "224.1",      "225.1", "226.1", "227.1", "228.1", "229.1", "230.1", "231.1",      "232.1", "233.1", "234.1", "235.1", "236.1", "237.1", "238.1",      "239.1", "240.1", "241.1", "242.1", "243.1", "244.1", "245.1",      "246.1", "247.1", "248.1", "249.1", "250.1", "251.1", "252.1",      "253.1", "254.1", "255.1", "256.1", "257.1", "258.1", "259.1",      "260.1", "261.1", "262.1", "263.1", "264.1", "265.1", "266.1",      "267.1", "268.1", "269.1", "270.1", "271.1", "272.1", "273.1",      "274.1", "275.1", "276.1", "277.1", "278.1", "279.1", "280.1",      "281.1", "282.1", "283.1", "284.1", "285.1", "286.1", "287.1",      "288.1", "289.1", "290.1", "291.1", "292.1", "293.1", "294.1",      "295.1", "296.1", "297.1", "298.1", "299.1", "300.1", "301.1",      "302.1", "303.1", "304.1", "305.1", "306.1", "307.1", "308.1",      "309.1", "310.1", "311.1", "312.1", "313.1", "314.1", "315.1",      "316.1", "317.1", "318.1", "319.1", "320.1", "321.1", "322.1",      "323.1", "324.1", "325.1", "326.1", "327.1", "328.1", "329.1",      "330.1", "331.1", "332.1", "333.1", "334.1", "335.1", "336.1",      "337.1", "338.1", "339.1", "340.1", "341.1", "342.1", "343.1",      "344.1", "345.1", "346.1", "347.1", "348.1", "349.1", "350.1",      "351.1", "352.1", "353.1", "354.1", "355.1", "356.1", "357.1",      "358.1", "359.1", "360.1", "361.1", "362.1", "363.1", "364.1",      "365.1", "366.1", "367.1", "368.1", "369.1", "370.1", "371.1",      "372.1", "373.1", "374.1", "375.1", "376.1", "377.1", "378.1",      "379.1", "380.1", "381.1", "382.1", "383.1", "384.1", "385.1",      "386.1", "387.1", "388.1", "389.1", "390.1", "391.1", "392.1",      "393.1", "394.1", "395.1", "396.1", "397.1", "398.1", "399.1",      "400.1", "401.1", "402.1", "403.1", "404.1", "405.1", "406.1",      "407.1", "408.1", "409.1", "410.1", "411.1", "412.1", "413.1",      "414.1", "415.1", "416.1", "417.1", "418.1", "419.1", "420.1",      "421.1", "422.1", "423.1", "424.1", "425.1", "426.1", "427.1",      "428.1", "429.1", "430.1", "431.1", "432.1", "433.1", "434.1",      "435.1", "436.1", "437.1", "438.1", "439.1", "440.1", "441.1",      "442.1", "443.1", "444.1", "445.1", "446.1", "447.1", "448.1",      "449.1", "450.1", "451.1", "452.1", "453.1", "454.1", "455.1",      "456.1", "457.1", "458.1", "459.1", "460.1", "461.1", "462.1",      "463.1", "464.1", "465.1", "466.1", "467.1", "468.1", "469.1",      "470.1", "471.1", "472.1", "473.1", "474.1", "475.1", "476.1",      "477.1", "478.1", "479.1", "480.1", "481.1", "482.1", "483.1",      "484.1", "485.1", "486.1", "487.1", "488.1", "489.1", "490.1",      "491.1", "492.1", "493.1", "494.1", "495.1", "496.1", "497.1",      "498.1", "499.1", "500.1", "501.1", "502.1", "503.1", "504.1",      "505.1", "506.1", "507.1", "508.1", "509.1", "510.1", "511.1",      "512.1", "513.1", "514.1", "515.1", "516.1", "517.1", "518.1",      "519.1", "520.1", "521.1", "522.1", "523.1", "524.1", "525.1",      "526.1", "527.1", "528.1", "529.1", "530.1", "531.1", "532.1",      "533.1", "534.1", "535.1", "536.1", "537.1", "538.1", "539.1",      "540.1", "541.1", "542.1", "543.1", "544.1", "545.1", "546.1",      "547.1", "548.1", "549.1", "550.1", "551.1", "552.1", "553.1",      "554.1", "555.1", "556.1", "557.1", "558.1", "559.1", "560.1",      "561.1", "562.1", "563.1", "564.1", "565.1", "566.1", "567.1",      "568.1", "569.1", "570.1", "571.1", "572.1", "573.1", "574.1",      "575.1", "576.1", "577.1", "578.1", "579.1", "580.1", "581.1",      "582.1", "583.1", "584.1", "585.1", "586.1", "587.1", "588.1",      "589.1", "590.1", "591.1", "592.1", "593.1", "594.1", "595.1",      "596.1", "597.1", "598.1", "599.1", "600.1", "601.1", "602.1",      "603.1", "604.1", "605.1", "606.1", "607.1", "608.1", "609.1",      "610.1", "611.1", "612.1", "613.1", "614.1", "615.1", "616.1",      "617.1", "618.1", "619.1", "620.1", "621.1", "622.1", "623.1",      "624.1", "625.1", "626.1", "627.1", "628.1", "629.1", "630.1",      "631.1", "632.1", "633.1", "634.1", "635.1", "636.1", "637.1",      "638.1", "639.1", "640.1", "641.1", "642.1", "643.1", "644.1",      "645.1", "646.1", "647.1", "648.1", "649.1", "650.1", "651.1",      "652.1", "653.1", "654.1", "655.1", "656.1", "657.1", "658.1",      "659.1", "660.1", "661.1", "662.1", "663.1", "664.1", "665.1",      "666.1", "667.1", "668.1", "669.1", "670.1", "671.1", "672.1",      "673.1", "674.1", "675.1", "676.1", "677.1", "678.1", "679.1",      "680.1", "681.1", "682.1", "683.1", "684.1", "685.1", "686.1",      "687.1", "688.1", "689.1", "690.1", "691.1", "692.1", "693.1",      "694.1", "695.1", "696.1", "697.1", "698.1", "699.1", "700.1"     ), class = "data.frame")


title: Which quantile type (quantile calculation method) is best for small samples?

text: The R help page of the quantile function lists 9 types of quantiles, i.e. 9 methods to calculate sample quantiles. Which method suits small sample sizes best? Which method will give the most similar quantile estimates if I take samples of sizes 6, 7, 8 of the same population? The distribution looks similar to log-normal; I'm most interested in percentiles 5, 50 and 95.


title: Using MLE to estimate the parameters of a strange exponential distribution

text: Is it possible to use MLE to estimate the parameters of a (very strange) exponential distribution? As you can see from the picture below, there is not much variation in the data. I wish to estimate the parameters that generate that form of exponential distribution![Exponential distribution](http://i.stack.imgur.com/dMvjH.jpg)


title: How to assess statistical significance of the accuracy of a classifier?

text: I have the classifier accuracy output in percentage and the number of input samples. Is there any test that can tell how statistically significant is the result based on this information. Thanks


title: How does one express the decrease in minimal type II error bound for each observation added?

text: **Problem** : I have a "classifier" that uses some arbitrary hypothesis test on observations from one of two known probability distributions:   1. $P_0$ (null hypothesis $H_0$) is a zero-mean Gaussian $\mathcal{N}(0,\sigma_0^2)$   2. $P_1$ (alternate hypothesis $H_1$) is a mixture of the zero-mean Gaussian from $P_0$ and some other zero-mean Gaussian with a higher variance $p\mathcal{N}(0,\sigma_0^2)+(1-p)\mathcal{N}(0,\sigma_1^2)$ where $\sigma_0^2<\sigma_1^2$ The parameters $\sigma_0^2$, $\sigma_1^2$, and $p$ are assumed to be known. Think of it as testing a machine. $H_0$ means that the machine is operating in good state, and $H_1$ means that the machine needs service. The "classifier" is allowed some fixed Type I Error $\alpha$. I am interested in determining how the _lower_ bound on Type II Error $\beta$ behaves as I use more observations. Specifically, I am interested in how the lower bound on $\beta$ decreases for each additional observation. Reason for this is that there is cost per each observation I make. **Question** : I am trying to find how much I can lower the _lower_ bound on $\beta$ with each additional observation, and knowing the parameters $\sigma_0^2$, $\sigma_1^2$, and $p$ of the distributions being tested. Can it be related to the sum of errors for a _single_ observation? **Prior knowledge** : For a single observation, I know that I can bound the sum of errors $\alpha+\beta\geq 1-TV(P_0,P_1)$ where $TV$ is total variation distance per question I asked at math.SE. The total variation between the two distributions in my problem is actually quite nice, as it's just $pTV(\mathcal{N}(0,\sigma_0^2),\mathcal{N}(0,\sigma_1^2))$. I also know that I can upper-bound $TV$ by the square root of half of Kullback-Leibler divergence using Pinsker's Inequality (there are additional benefits to using KL bound which are out of scope of this question...) However, I do not know how to relate the bound on one observation to lowering of the bound with an additional observation. **My mathematical sophistication** : I have graduate-level background in probability theory (the M.S. Engineering course). I also have graduate-level background in information theory and signal theory. Unfortunately, I have very limited knowledge of measure theory...


title: Covariate modeling with a within-subject factor

text: A researcher collected some data of a response/outcome variable Y from, for example, 20 subjects. In addition to an explanatory (within-subject or repeated-measures) variable A (a factor with two levels, A1 and A2), he would also want to consider another explanatory (and continuous) variable X. He's interested in finding out the effect of factor A on Y, and in the meantime would like to control for the variability in Y due to X. However, the subtle issue here is that X is correlated with Y to some extent, and, more importantly, the two levels (A1 and A2) may have different average value of X. The complication of the latter fact in modeling is that the effect of X on Y may partially or even fully explains the effect of A. I know this may be a controversial issue because X is correlated with the levels of A to some extent, but the investigator wants to account for the variability **within** each level of A. If I construct a model like               Y ~ X + A      or,               Y ~ X*A      with lme4 of nlme package in R, I believe that the result of A effect from anova() should be interpreted with X taking the overall mean value, but that would not allow me to obtain the effect of A while 'controlling' for X (considering the fact that X has different mean value across the two levels of A; i.e., the average value of X at A1 is different from the average value of X at A2)! I've thought of scaling Y by X with a model like               Y/X ~ A      but I don't feel so comfortable with it. Any suggestions on modeling this situation? Is there such a thing as multivariate LME?


title: How to digest statistical context?

text: Firstly, I suppose that **not** all active members of this interesting site are statisticians as their job. Otherwise the question being asked as follows does not make any sense! I respect them of course. but I need an explanation that is a bit more practical rather than conceptual.   I start with an example from Wikipedia to define `point process`: > Let S be locally compact second countable Hausdorff space equipped with its > Borel σ-algebra B(S). Write $\mathfrak{N}$ for the set of locally finite > counting measures on S and $\mathcal{N}$ for the smallest σ-algebra on > $\mathfrak{N}$ that renders all the point counts ... measurable. To me this has no meaning. **An explanation in an engineering context is more understandable to me.**   **Comment:** Most of time I found Wikipedia's explanations **useless** due to similar **complicated text** (at least to me). From my experience there are only two types of reference books for statistics: **a)** extremely simplified **b)** extremely complicated!   Reading both has no benefits for me at all! ### Question:   * Do you have a solution for this problem? or similar experience ...? ### Important Note: Please answer if you really understood what the question is about. Please don't recommend simple solutions that could be easily found via googling! None of them is of my answer. ** For those who found this post useful there are benefits to check also this which discusses a related topic from different prospective.


title: Data prep / variable creation for predictive models

text: I was reading a couple of the write ups from a Kaggle challenge: Here is one and another and it got me wondering about variable creation in data mining and why there seems to exist so few texts or case studies available. Are there any example out there besides Dorian Pyle book? I am thinking full case studies of prediction challenges and feature creation / variable deriving and how one should approach it: domain knowledge of course, but also from EDA.


title: Estimator of mu for gaussian variable with penalty

text: I have a question about the following problem: estimate $\mu=\mu_1,...\mu_n$ when $Y_i \sim N(\mu_i,\sigma^2)$ using a ridge like penalty $$ \min_\mu \sum_i(Y_i-\mu_i)^2 + \lambda\sum_i\mu_i^2 $$ I have calculated the closed form decomposition of the generalized cross validation error and I know the optimal $\lambda$ depends on $\mu$ and $ \sigma$. My question is how calculate the optimal $\lambda$ and the cross validation error with a cross validation procedure ? I don't see how to do that. Thanks for your help


title: Manipulating ordinal variables to segment data

text: Im working with a negative binomial regression where the dependent variable is number of trips by any given mode (i.e. car trips, train trips etc). I would like to create a variable to segment the data based on theoretical grounds. The variables that I would like to transform to do so are a set of six questions regarding individual preference and measured on a five point likert scale. I would like to segment the data into two mutually exclusive groups according to two different types of residential preference (let`s say Group A and Group B). In this case, and again, based on the field`s literature, three of the questions relate to Preference A and three questions relate to Preference B. Now, since I do not have enough factors to conduct factor analysis, I was thinking of coming up with a rule to decide whether an individual fits in one of the two classes. My thinking is as follows: 1.Take the median of the first group of three questions (related to Preference A) for the ith individual. 2.Take the median of the second group of three questions (related to Preference B) for the ith individual. *Since the mean cannot be used as a centrality measure for ordinal data because of the unknown distance between each point, im assuming that even if the distances might be different for two different individuals, for any given individual i the distance from 1 to 2 and from 3 to 4 are the same for ALL six questions.  3.In that sense, if the median of the Preference A questions is higher than the median for the Preference B questions, I am assigning that indiviudal to preference group A and viceversa. In the case the median is the same for both sides I am assuming that the individual has no particular preference for either group. Of course I undesrtand it is a very strong assumption given the we do not know the distance between each point but given the constraints of a secondary dataset, I am wondering to what extent this reasoning is adequate to address the problem. Thanks in advance.


title: How likely am I to be descended from a particular person born in the year 1300?

text: In other words, based on the following, what is p? In order to make this a math problem rather than anthropology or social science, and to simplify the problem, assume that mates are selected with equal probability across the population, except that siblings and first cousins never mate, and mates are always selected from the same generation.   * $n_1$ -- initial population   * $g$ -- the number generations.   * $c$ -- the average number of children per couple. (If necessary for the answer, assume that every couple has exactly the same number of children.)   * $z$ -- the percentage of people who have no children, and who are not considered part of a couple.   * $n_2$ -- population at the final generation. (Either $n_2$ or $z$ should be given, and (I think) the other can be calculated.)   * $p$ -- probability of someone in the final generation being a descendant of a particular person in the initial generation. These variables can be changed, omitted, or added to, of course. I am assuming for simplicity that $c$ and $z$ do not change over time. I realize this will get a _very_ rough estimate, but it's a starting point. **Part 2 (suggestion for further research):** How can you consider that mates are not selected with globally uniform probability? In reality, mates are more likely to be of the same geographical area, socio-economic background, race, and religious background. Without researching the actual probabilities for this, how would variables for these factors come into play? How important would this be?


title: Weighted average semi-parametric regression in R

text: Are there any packages that supports weighted average semi-parametric regressions in R? An example of such a regression is in the links below. I see that there is package GAM in R for Generalized Additive Models. However, it is not clear to me whether this procedure is the same as weighted average semi-parametric regression. There is also a package called SPM Some background on semi-parametric methods as used by Connor is here. This paper refers extensively to improvements in the methodology vs. a prior paper by Connor here Also, attached is a paper that discusses non-parametric regressions in R Thank you


title: Effect size calculations in mixed ANOVA design for power analysis

text: I need to calculate effect size for a mixed ANOVA design. I need it because I'm using GPower3 to calculate a rough sample size for a full study. I have pilot data. Based on this page: http://www.xlstat.com/en/features/statistical-power-analysis-of-variance- anova.htm and the GPower3 software, it looks like there are 3 ways to calculate an effect size:   * eta²   * use variances   * use group means Is there any rational for picking one over another? Since I have the pilot data, I guess I was just going to use the variances to calculate the effect size. The formula I have for that is: f (effect size) = sqrt(Variance Explained / Variance Error) Based on the summation rule: Total Variance = Variance Explained + Variance Error Is this correct? If I have unequal group sizes, could I use the pooled variance formula for the Variance Error term? I also noticed that this effect size formula using the variances is quite similar to the F statistic: F = Variance Explained / Variance Error Is this correct and/or is there a relationship between them? Thanks in advance.


title: The larger $R^2$ the better?

text: I want to show that the variable $X$ is significant.   In model 1, I use financial statements variables, other macroeconomic variables and $X$.   Here, $X$ is siginificant at a 10% level and $R^2$ is 0.731.   Since financial statements could have time lag, I have to check that effect so I made adjusted model 2.   In here $X$ is insignificant but $R^2$ is 0.721.   So, I can say that "since $R^2$ is decreased, model 1 is better and $X$ is significant variable"?


title: Finding all largest sequences

text: What would be a good / efficient algorithm or approach to find all the largest sequences within a list of chains with varying lengths? For instance these chains:               0: [4,5,2,1,3,6,8,9]     1: [2,1,3,4,5,6,7,8,10,11]     2: [15,12]      The longest common subsequence would be [2,1,3,6,8] but there is another important subsequence [4,5,6,8] that we would be needed to find. All other sequences with 2 or more elements would be inside those two chains. For instance [1,3,6] or [5,6,8] or [2,1] are all inside the two resulting chains. Note: the example uses integers for more clarity but the real data will use floats, and they will be taken as a "match" when both are separate within a given threshold.


title: Discrepancy between coefficient and mean difference in predicted values of logistic regression

text: I'm using a poisson-binominal logsitic regression model to analyze a list experiment (item count technique) where the outcome variable is a binary response of the respondent to a sensitive item (e.g., whether or not a black family moving next door to you will make you angry). I'm using a poisson- binomial distribution rather than the binominal distribution because the probability of giving an affirmative answer to control items differs across items, hence I'm looking at the sum of independent but heterogenous bernoulli random variables, each with different success probabilities. I'm interested in explanatory variable A (also dichotomous). From this model, the estimated coefficient of explanatory variable A is -2.275 and the standard error is 0.972 (i.e., the coefficient is significant at the 95% level). However, when I calculate the mean difference in predicted values across two datasets where explanatory variable A = 1 and A= 0 (all other covariates are set to their observed values), the mean difference is -0.093 and the standard error is 0.287 (i.e., the 95% confidence intervals cross 0). What would explain this difference between a significant coefficient estimate and a not significant mean difference in predicted values?


title: Dealing with linear dependent variables

text: I have a large dataset with many subject each with responses from a consecutive year going back 10 years (ie 100,000 persons per year (not necessarily 10 data points per person as they may not have been part of the study in prior years) dating back 10 years). I have data on each specific person that varies within each year, and data that specifically varies with just each year, meaning that each of the 100,000 person in the same year will have the same value if they belong to a certain year(and that value changes with each year). Essentially I have just 10 values in that variable. What complicates this is I have multiple variables like this (housing values, opinion surveys, etc) that have a significant impact on how that individual behaved given that year. Also I have made a few lag variables for each of those variables that is a change in value from a specific time in the past (ie change in home value from 1 yr ago, 2 years ago, 4 years ago, etc). I am running into a linear dependence issue with some of these variables, and when I run the regression, the model statistics are the same, but I can get very different results as well as predictions on test data. What is the best way to handle these types of variables by transformations, modeling techniques, or other methods?


title: Clustering on a data set with mixed variables

text: I have a data set consisting of $n$ elements with $d$ features for each element ($x_{i,f}$ means the value for the f-th feature of the i-th element). I would like to cluster this data set into $k$ clusters. One problem I have is that one feature is a nominal one and another is a discrete ordinal one. As an example consider elements which have the following features:   * $x_{i,1}$ : the height of a person   * $x_{i,2}$ : the weight of a person   * $x_{i,3}$ : country where the person lives   * $x_{i,4}$ : nr. of friends the person has Is it ok to use a simple k-means algorithm with an euclidian distance measure? I would introduce a indicator variable $\delta$ with the following meaning for feature $x_{i,3}$: $$\delta(i,j) = \begin{cases} 0, & \text{if }x_{i,3} = y_{j,3}\text{,}\\\ 1, & \text{else.}\end{cases}$$ So two objects $i$ and $j$ have no distance for their third feature if this feature is the same (same country) and 1 otherwise. Or do you know a better way to do a cluster analysis in this case?


title: Which statistical test to use? I am really confused.

text: I have a research paper; I used two scales to measure marital adjustment and social roles. I applied them to 150 people. I would like say something about how males who got high scores on the social roles scale compare to women who also got high scores on marital adjustment scores, or vise versa. The scores are both interval measures. I'm confused what tests to use. I want to test a hypothesis that if men get higher social role scores and women get higher on social role scores, then marital adjustment scores are also higher.


title: ARMA modeling in R

text: I have a time series. I want to model it using ARMA, which will be used for forcasting. In R I am using `arima()` function to get the coefficients. But `arima()` requires order(p,d,q) as input. What is the simplest way in R to arrive at a good value for p and q (with d = 0) so that I don't overfit?


title: Understanding shape and calculation of confidence bands in linear regression

text: I am trying to understand the origin of the curved shaped of confidence bands associated with an OLS linear regression and how it relates to the confidence intervals of the regression parameters (slope and intercept), for example (using R):               require(visreg)     fit <- lm(Ozone ~ Solar.R,data=airquality)     visreg(fit)      ![enter image description here](http://i.stack.imgur.com/Q7mQK.png) It appears that the band is related to the limits of the lines calculated with the 2.5% intercept, and the 97.5% slope, as well as with the 97.5% intercept, and the 2.5% slope (although not quite):               xnew <- seq(0,400)     int <- confint(fit)     lines(xnew, (int[1,2]+int[2,1]*xnew))     lines(xnew, (int[1,1]+int[2,2]*xnew))      ![enter image description here](http://i.stack.imgur.com/2KLTH.png) What I don't understand are two things:   1. What about the combination of 2.5% slope & 2.5% intercept as well as 97.5% slope and 97.5% intercept? These give lines that are clearly outside the band plotted above. Maybe I don't understand the meaning of a confidence interval, but if in 95% of the cases my estimates are within the confidence interval, these seem like a possible outcome?   2. What determines the minimum distance between the upper and lower limit (i.e. close to the point where the two lines added above intercept)?  I guess both questions arise because I don't know/understand how these bands are actually calculated. How can I calculate the upper and lower limits using the confidence intervals of the regression parameters (without relying on predict() or a similar function, i.e. by hand)? I tried to decipher the predict.lm function in R, but the coding is beyond me. I'd appreciate any pointers towards relevant literature or explanations suitable for stats beginners. Thanks.


title: Generalized linear model with mean contrast pairwise comparison versus chi-square tests

text: I have binary data on prevalence of infection (infected/all analysed) and I want to compare this proportion between groups (species and sexes). I would like to ask if by conducting a Generalized Linear Model with binomial distribution and mean contrast pairwise comparison (with 1/0 as target variable and the groups and their interaction as factors), the pairwise results I get are "the same" as conducting separate chi-Square tests. If it is the same, I am using SPSS and I get a significance for the model effects and a slightly different significance for the pairwise comparisons. Which one should I report? What do these differences mean? If it is not the same, I would like to ask what is the correct way to construct the contingency table for comparing, for example, species A with 12/30 with species B with 8/22. And within species A, males have 10/17 and females 2/13, while within species B, males have 5/10 and females with 3/12. Can I use chi-square with different sample sizes of each group? Does anyone know how to do this in SPSS (I think it requires same sample size)? For Species A versus Species B, the result from chi-square is Pearson=0.07, P=0.79. Thus no differences between species. For Male_A versus Female_A, result from chi-square is Pearson=5.79, P=0.02. Thus differences between sexes. Is this correct? Now imagine for A I have 4/29 and for B I have 2/16. If I am doing the chi- square table correctly, I have an additional problem of less than 5 in some cells, for which I only get Fisher's exact test. The results are only for Fisher Exact Probability Test, P one-tailed=0.64 and P two-tailed=1. Does this mean there are no differences? Which P is relevant and how can I get a chi- square value in these situations?


title: Universal Approximation Theorem — Neural Networks

text: I have posted this question elsewhere--MSE-Meta, MSE, TCS, MetaOptimize. Previously, no one had given a solution. But now, here is a really excellent and comprehensive answer. Universal approximation theorem states that "the standard multilayer feed- forward network with a single hidden layer, which contains finite number of hidden neurons, is a universal approximator among continuous functions on compact subsets of Rn, under mild assumptions on the activation function." I understand what this means, but the relevant papers are too far over my level of math understanding to grasp why it is true or how a hidden layer approximates non-linear functions. So, in terms little more advanced than basic calculus and linear algebra, how does a feed-forward network with one hidden layer approximate continuous, non- linear functions? The answer need not necessarily be totally concrete.


title: What, if any, dissimilarity is preserved in partial least squares (PLS)?

text: When we perform a principal components analysis (PCA) on a multivariate data set we are interested in finding orthogonal components that explain maximal variance in the data set. We can form a biplot of the data using the scores and the loadings, and the locations of the sample points in the biplot are an approximation of the Euclidean distance between the samples. In PLS, we extract orthogonal components from a predictor data set that have maximal covariance with the response (vector or matrix). We also get scores and loadings as part of the analysis and can draw a biplot of these scores. What, if any, dissimilarity is represented by the Euclidean distances _on the biplot_ between sample points? One of the reasons I ask is that with PCA, we can apply a transformation to the data prior to applying PCA such that the Euclidean distance between samples on the biplot approximates the Euclidean distance between samples _in_ the transformed data, but in the untransformed data the distance represented is some other distance. For example, by applying the Hellinger transformation (rows are standardised by their row sum and then a square root transformation is applied) to the raw data, a PCA applied to the transformed data will reflect the Hellinger distances between the observations. I wonder if a similar principal might hold for PLS?


title: How to estimate parameters for a Kalman filter

text: In a previous question, I inquired about fitting distributions to some non- Gaussian empirical data. It was suggested to me offline, that I might try the assumption that the data is Gaussian and fit a Kalman filter first. Then, depending on the errors, decide if it is worth developing something fancier. That makes sense. So, with a nice set of time series data, I need to estimate **several** variable for a Kalman filter to run. (Sure, there is probably an R package somewhere, but I want to actually **learn** how to do this myself.)


title: Quantifying how successful results of clustering are with a test

text: I have a three types of experiments with many measurements taken. I want to show that using clustering on the unlabelled measurements, I can uncover the three groups. Now using kmeans plus, I get around 60% accuracy. Which if it was binary, would be poor, but how do I quantify how good that is for three categories? Is there a p value I can compute and compare against? I am thinking that if there were 100 categories 60% for the correct labelling is great and far from random, so how can I do this with a test?


title: What to do when Kolmogorov-Smirnov test is significant for residuals of parametric test but skewness and kurtosis look normal?

text: I have conducted a parametric test in a study, n=290. I want to assess whether the residuals of this test are normally distributed.   * The skewness and kurtosis of the residuals are -0.017 and -0.438 respectively. I think this is considered as normal.    * Unfortunately, the Kolmogorov-Smirnov of the residuals has a p-value of 0.021. I think this is considered not normal.  ### Question **What should I do when skewness and kurtosis look normal but Kolmogorov- Smirnov is significant?**


title: Model Building: Missing Data or Large Gap between data points

text: I am currently trying to build a model using a data set that has large gap between data points. When I look for the correlation I clearly see a negative regression line. But I am worried about the gap that exist between the poins. I build a simple linear model though this has high R squared I don't think simple linear regression is the best model to that fits the data. This looks like it has a negative exponential behavior. I thought to post here to get some expert thoughts on what I should do when you deal with the data that has a large gap between points and does this data has a linear relationship or strong non linear relationship? Data Set:                  density  co2     1     20.4 38.8     2     27.4 31.5     3    106.2 10.6     4     80.4 16.1     5    141.3  7.7     6    130.9  8.3     7    121.7  8.5     8    106.5 11.1     9    130.5  8.6     10   101.1 11.1     11   123.9  9.8     12   144.2  7.8     13    29.5 31.8     14    30.8 31.6     15    26.5 34.0     16    35.7 28.9     17    30.0 28.8     18   106.2 10.5     19    97.0 12.3     20    90.1 13.2     21   106.7 11.4     22    99.3 11.2     23   107.2 10.3     24   109.1 11.4      Plot: ![Simple Linear Regression Plot](http://i.stack.imgur.com/MWjgD.png) Summary of Linear Model:               Coefficients:                 Estimate Std. Error t value Pr(>|t|)         (Intercept) 38.12948    1.21768   31.31  < 2e-16 ***     density     -0.24247    0.01261  -19.22 3.04e-15 ***      In addition if I transfer both density and co2 as log transform variables, then I see following behavior. Since data is missing at the middle its really hard to stick to a log transformed model or the base model. ![enter image description here](http://i.stack.imgur.com/pPUZN.png)


title: Are there any libraries available for CART-like methods using sparse predictors & responses?

text: I'm working with some large data sets using the gbm package in R. Both my predictor matrix and my response vector are pretty sparse (i.e. most entries are zero). I was hoping to build decision trees using an algorithm that takes advantage of this sparseness, as was done here). In that paper, as in my situation, most items have only a few of the many possible features, so they were able to avoid a lot of wasted computation by assuming that their items lacked a given feature unless the data explicitly said otherwise. My hope is that I could get a similar speedup by using this sort of algorithm (and then wrapping a boosting algorithm around it to improve my predictive accuracy). Since they didn't seem to publish their code, I was wondering if there were any open-source packages or libraries (in any language) that are optimized for this case. Ideally, I'd like something that could take a sparse matrix directly from R's `Matrix` package, but I'll take what I can get. I've looked around and it seems like this sort of thing should be out there:   * Chemists seem to run into this issue a lot (the paper I linked above was about learning to find new drug compounds), but the implementations I could find were either proprietary or highly specialized for chemical analysis. It's possible one of them could be re-purposed, though.   * Document classification also seems to be an area where learning from sparse feature spaces is useful (most documents don't contain most words). For instance, there's an oblique reference to a sparse implementation of C4.5 (a CART-like algorithm) in this paper, but no code.   * According to the mailing list, WEKA can accept sparse data, but unlike the method in the paper I linked above, WEKA isn't optimized to actually take advantage of it in terms of avoiding wasted CPU cycles. Thanks in advance!


title: GLMM predictor variables

text: I have a question regarding predictor variables in a GLMM analysis. I am running a study investigating multi-modal communication in primates, specifically, primate gestures that accompany vocalizations. I have measured the call rate, call duration and peak frequency (response variables) for 200 calls (100 accompanied by gestures, 100 without gestures). I would like to use the following predictor variables: 1) Call type, 2) Presence of gesture, 3) Gesture type. Is it possible to use GESTURE TYPE here as a predictor variable? I am not sure if it's appropriate since this will only be for 50% of the overall the data (where calls are accompanied by gestures).


title: Borrowing strength

text: What are the principles of Borrowing Strength? What does it mean in terms of estimating parameters for hierarchical models? Where can this information can be read from?


title: logistic regression with dummy variables for fractional factorial design

text: We have conducted a survey experiment with varying amounts of incentive (factor 1 = I1, I2, I3, I4, I5). The experiment was conducted stepwise in three subsequent studies (factor 2 = S1, S2, S3). Unfortunately, the factors are not fully crossed. The experimental design is as follows: ![enter image description here](http://i.stack.imgur.com/WPUjE.jpg) Outcome variable of the planned analysis is recruitment for a panel survey (y=0,1) dependent on incentive group. I was wondering if it is statistically correct to include all five incentive groups and furthermore controlling for the point in time in one regression analysis. The five incentive groups should be included as dummy variables in the regression analysis. When constructing interaction terms of incentive group and study, one is confronted with the problem of multicollinearity. E.g., respondents of study 1 (S1) have no chance to be in incentive groups 3 to 5 by design. In the table below I have illustrated the data matrix with all 8 possible patterns. The model can be formally described as follows: ![enter image description here](http://i.stack.imgur.com/gvEa3.jpg) ![enter image description here](http://i.stack.imgur.com/hgKBG.jpg) Is it correct to leave out incentive group 1 (I1) and study 1 (S1) as I have done? Finally, I was wondering if the varying numbers of respondents per experimental group will be an issue besides less precise estimates for the groups with smaller numbers. Thank you very much for potential answers!


title: Program yourself or use data mining toolkit?

text: I'm quite comfortable using Python and additional ML/Visualization packages for ML and data mining. How do the free data mining software RapidMinder, KNIME, etc. compare to this approach, and do they have a practical advantage over using ML libraries? Are these big software packages more than fancy drag- and-drop which slows you down if you have programming skills? What is the approach of people who do serious data mining and who possess programming skills to choose? Do they prefer to write code?


title: Probit model marginal effects

text: p<-glm(GRADE~GPA+TUCE+PSI,family="binomial"(link="probit")); summary(p)      is the probit model               Coefficients:                   Estimate    Std. Error  z value   Pr(>|z|)        (Intercept)   -7.45231    2.57152     -2.898    0.00376 **     GPA            1.62581    0.68973      2.357    0.01841 *      TUCE           0.05173    0.08119      0.637    0.52406       PSI            1.42633    0.58695      2.430    0.01510 *       If we have mean of TUCE and PSI=1, what is the marginal effect of GPA on Pr(GRADE=1)?               pnorm(-7.45231+1.62581*mean(GPA)+0.05173*mean(TUCE)+1.42633*1)*1.62581       dnorm(-7.45231+1.62581*mean(GPA)+0.05173*mean(TUCE)+1.42633*1)*1.62581        I'm wondering which one of `pnorm` and `dnorm` is correct to use if you want the marginal effect.


title: test the significance of clusters

text: Good morning, I am analyzing a dataset composed by 364 subjects and 13 binary variables (0,1 = absence,presence). I am testing possible association (co- presence) of my variables. To do this, I was trying with cluster analysis. My main interest is to check for the significance of the obtained clusters. First, I tried with the `pvclust()` function, by using `method.hclust="ward"` and `method.dist="binary"`. Altoghether it works (clusters and significance obtained). However, I'm not convinced by the distance matrix. Association between variables are indeed different from results obtained in PAST by using Ward on a Jaccard matrix (that should be ok for binary data). Moreover, when I try to obtain a Jaccard matrix in R from my data, by using the Vegan package               mydistance<-vegdist(t(data),method="jaccard")      I receive the following error message:               Error in rowSums(x, na.rm = TRUE) : 'x' must be numeric      below an subset from my dataset:               variable1   variable2   variable3   variable4   variable5   variable6   variable7   variable8   variable9   variable10  variable11  variable12  variable13     case1   0   0   0   0   0   1   0   0   1   1   0   0   0     case2   0   0   0   0   0   1   0   NA  NA  1   0   0   0     case3   0   0   0   0   0   1   0   0   1   1   0   0   0     case4   1   0   0   0   0   1   0   1   0   1   0   0   0     case5   0   0   0   0   0   1   0   0   1   1   0   0   0     case6   0   1   0   0   0   1   0   1   0   1   0   0   0     case7   0   1   0   0   0   1   0   0   1   1   0   0   0     case8   0   0   0   0   0   1   0   1   0   1   0   0   0     case9   0   0   0   0   0   1   0   1   0   1   0   0   0     case10  0   0   0   0   0   1   0   0   1   1   0   0   0     case11  1   0   0   1   0   1   1   1   0   1   0   0   0     case12  0   0   0   1   1   0   1   1   0   1   0   0   0     .....         So, my questions are the following: Is the Jaccard index a good strategy for my kind of data? Is binary distance used in pvclust is theoretically more correct? Is there any alternative to pvclust for testing the significance of my clusters? Thanks in advance Marco


title: How to tell quantitatively whether 1D data is clustered around 1 or 3 values?

text: I've got some data on the time between heart beats of a human. One indication of ectopic (extra) beats is that these intervals are clustered around three values instead of one. How can I obtain a quantitative measure of this? I'm looking to compare multiple data sets, and these two 100-bin histograms are representative of all of them. ![enter image description here](http://i.stack.imgur.com/fyPFU.png) I could compare the variances, but I want my algorithm to be able to detect whether there is one or three clusters in each case without comparing to the other cases. This is for offline processing, so there's a lot of computation power available, if that's needed.


title: Variance of the Kaplan-Meier estimate for dependent observations

text: Can someone help me find a way to estimate the variance of the Kaplan-Meier estimate with dependent observations? Specifically, I have failure time data from patients with several different observations for each patient (and different patients may have different number of observations). The observations for different patients are assumed to be independent but observations from the same patient are expected to be dependent. I was suggested two publications, "Kaplan-Meier Analysis of Dental Implant Survival: A Strategy for Estimating Survival with Clustered Observations" and "The Kaplan-Meier Estimate for Dependent Failure Time Observations", the second of which is cited by the first. However I was unable to make sense of these. The first has errors and the second seems far more rigorous but the equation for the variance does not make sense (double integration on a 1-form).


title: How to deal with label switching issues with Matlab's classification trees?

text: I have a situation that seems trivial but I can't figure it out. I have a dataset in Matlab that has categorical values. For example:               Outlook,Temperature,Humidity,Windy,Play     sunny,hot,high,false,no     sunny,hot,high,true,no     overcast,hot,high,false,yes     rainy,mild,high,false,yes     rainy,cool,normal,false,yes     rainy,cool,normal,true,no     overcast,cool,normal,true,yes      In order to use `ClassificationTree` (or other algorithm in machine learning), I need to convert the dataset to a matrix. I use (assuming that data contains the dataset):               double(data)      which assigns numbers to the values. Let's say `rainy=1`, `overcast=2`, `sunny=3` and so on. I train the ClassificationTree and it works fine. But here is my problem. If I want to predict on test data:               overcast,cool,normal,false,yes      I need to know the numbers that `double()` assigned to each of the categories. i.e. `overcast=2`. Using `double()` on the test set does not work because there is no guarantee that the numbers assigned by `double()` are the same to those assigned to the training dataset. I have found a really twisted way of doing it. For example, for `Outlook`:               d = zeros(size(test));     ls = getlevels(training.Outlook);     n = size(ls,2);     for i = 1:n         d(test.Outlook == ls(i),1) = i;     end      This will assign the numbers correctly because I check against each of the values for `Outlook` that I found in the training dataset. This is far from elegant, and there has to be something better. Any help would be appreciated.


title: Tobit Model with a Variable Lower Censoring Bound?

text: I am working with NBA contract data and about 20 percent of my free agent sample do not receive contracts. Additionally, the distribution for those who do receive contracts is exponential. For each individual who did not receive a contract, I know what the league-mandated minimum salary would have been. It is a function of Years of Service, and thus I know the assessed contract value of the player is between 0 and the mandated minimum. I am trying to estimate a regression equation on the full sample, and would like to use a tobit specification as my data is clearly left-censored. The problem is that the censoring bound is a variable that is a function of the Years of Service. Does anyone know of examples of a tobit model with a variable bound? Even better would be code to run such a model in Stata, but I'm not getting my hopes up.


title: How to find the inverse of non-standard CDF in Matlab?

text: The cumulative distribution function (CDF) for a truncated log-normal distribution is represented as shown below. How can I find its inverse in Matlab? `CDF_truncated_lognormal_distribution = lognormal_cdf(r)./(lognormal_cdf(rmax) - lognormal_cdf(rmin))` Where, `lognormal_cdf()` in above equation is a CDF of the standard log-normal distribution, and `rmax` and `rmin` are constants. I need to find the inverse of above CDF in order to plug uniform random numbers (=U[0,1]) as a CDF value in the inverse equation and determine random number `r` of the truncated log-normal distribution.