title: How can I pool bootstrapped p-values across multiply imputed data sets?

text: I am concerned with the problem that I would like to bootstrap the p-value for an estimate of $\theta$ from multiply imputed (MI) data, but that it is unclear to me how to combine the p-values across MI sets. For MI data sets, the standard approach to get to the total variance of estimates uses Rubin's rules. See here for a review of pooling MI data sets. The square root of the total variance serves as a standard error estimate of $\theta$. However, for some estimators the total variance has no known closed form or the sampling distribution is not normal. The statistic ${\theta}/{se(\theta)}$ may then not be t-distributed, not even asymptotically. Therefore, in the complete data case, one alternative option is to bootstrap the statistic to find variance, a p-value and a confidence interval, even if the samling distribution is not normal and its closed form unknown. In the MI case there are then two options:   * Pool the bootstrapped variance across MI data sets   * Pool the p-value or confidence bounds across MI data sets The first option would then again use Rubin's rules. However, I believe that this is problematic, if $\theta$ has a non-normal sampling distribution. In this situation (or more generally, in all situations) the bootstrapped p-value can be used directly. However, in the MI case, this would lead to multiple p-values or confidence intervals, which need to be pooled across MI data sets. So my question is: how should I pool multiple bootstrapped p-values (or confidence intervals) across multiply imputed data sets? I would welcome any suggestions on how to proceed, thank you.


title: Factor analysis problems

text: I am completing a factor analysis with 105 items and am having some problems. When I look at the scree plot for the "elbow" it suggests that there a 5 factors, however these five factors only account for 43% of all the variance. Looking at the eigen values it suggests that there are 23 factors…. The items are likert scales, with the respondents rating themselves as to how they react to given situations. It was planned that this would allow for the formation of a value that would indicate a propensity to respond poorly to situations, made up of several factors. The 105 items appear to be reliably assessing this (crombach's alpha is over .85). Is there any way of improving the analysis? When I look at the factor loadings of the items I can only really allocate the items to perhaps 3-4 factors (with one very large factor and multiple items that don't clearly describe any factors). Is there a way of weeding problematic items out and then redoing the analysis?


title: Testing for autocorrelation in standardised squared residuals of GARCH models

text: I was wondering at what lags should I perform the Ljung Box test on on standardised squared residuals. I used fGarch and the tests were done at lag 5,10 and 15 respectively with no evidence to reject the null hypothesis. However, I tried using rugarch as it has functions that fGarch doesn't have but the ljung box test were done at lag 1,3 and 7 with the null hypothesis being rejected. I'm getting confused as to which lags should be tested as it seems that different lags yields different response. Also isn't there a way to tell to rugarch which lag should be tested?


title: Binomial hierarchical model

text: I want to analyse the following data: 12 populations have been followed and a certain illness in the population is measured. I am interested in the probability of being ill. Now I divide these populations in males and females which would have their own probabilities in turn (I could in theory slice and dice the data differently). I would have intuitively expected that if a population is small its probability should be closer to the total population probability and vice versa. So this is why I believe a hierarchical model might be useful. To formalise, I have the following model in mind: > $y_{ij}$ ~ Binomial($N_{ij}$, $p_j$) for i=1..12, j=1 (total population), 2 > (male) and 3 (female) > > $p_1$ ~ Beta(1,1) > > $logit(p_j)$ ~ $logit(p_1)$ + $\epsilon$ for j=2,3 > > $\epsilon$ ~ Normal(0,$\sigma^2$) I wonder if this specification makes sense. Could you please point to some literature or other more suitable specifications?


title: k-means vs k-median?

text: I know there is k-means clustering algorithm and k-median. One that uses the mean as the center of the cluster and the other uses the median. My question is: when/where to use which?


title: Intuition behind the Stein's paradox

text: I had read wiki and some sources. jmanton's blog, Wasserman's blog * * * the background is that: You have Xi ∼ N(θi, 1), and we want to estimate the each θi.   Where Xi are `independent`。 For `MSE` risk of vector (X1,...,Xn), when n>=3.   * the MLE estimator $\hat θ_i = x_i$ is inadmissable .    * While James-Stein estimator is better than MLE in such case. * * * ## Looking from Shrinkage quote from Larry Wasserman: > Note that the James-Stein estimator `shrinks` {X} towards the origin. (In > fact, you can shrink towards any point; there is nothing special about the > origin.) I know the James-Stein estimator is special to MLE for its shrinkage behavior. I still dont get it, the variables are independent. **Why James-Stein estimator`shrinkage` to arbitrary point can improve the MSE risk versus where MLE do `not shrinkage` at all ??** * * * and seeing this from Larry Wasserman: > This can be viewed as an empirical Bayes estimator .... If view it from Empirical Bayes case, it is shrinkage to a overall mean, but that is nonsense if you have independent variable to estimate. * * * Is there better example or intuitive explanation for it?? THANKS.


title: A regression model whose response variable is the day of year that an annual event (usually) occurs

text: In this particular case I'm referring to the day on which a lake freezes. This "ice-on" date only occurs once a year, but sometimes it doesn't occur at all (if the winter is warm). So on one year the lake may freeze on day 20 (january 20th), and another year it might not freeze at all. The goal is to figure out drivers of ice-on date. Predictors would be things like fall/ winter air temperature each year. Year could be a predictor for the long-term linear trend. 1) Is the integer "day of year" a reasonable response variable (if not, what is?)? 2) How should one handle the years when the lake never froze?


title: Convergence of Gini index

text: > Let $\theta(F)=2\int^1_0(t-q_F(t))dt$, where $\displaystyle > q_F(t)=\frac{\int_0^tF^{-1}(s)ds}{\int_0^1F^{-1}(s)ds}$. For discrete distributions I'm assuming that $F^{-1}(s)=\inf\\{x:F(x)=s\\}$ (the quantile function). The exercise next defines a discrete distribution function in a strange way. > Let $F_p$ have probability $p$ at $x=0$ and $1-p$ at $x>0$. (ipsis verbis from the book) This doesn't seem to be a proper distribution, no? Then they ask me to > show $\theta(F_p) \xrightarrow[]{p\rightarrow 0}0$ and $\theta(F_p) > \xrightarrow[]{p\rightarrow 1}1 $. My problem begins immediately at finding $F^{-1}$. It seems to be always zero...


title: Favorable properties for penalized likelihood estimators

text: I'm reading through Fan and Li's SCAD paper on Nonconcave Penalized Likelihood. The idea is to penalize the loss function by $p_\lambda(\theta)$ where $\lambda$ can be some tuning parameter. Fan proposes that a proper penalty function should produce estimators with the following 3 properties, and also lists sufficient conditions for the penalty function to have these properties: * * * Unbiasedness: The resulting estimator is nearly unbiased when the true unknown parameter is large to avoid unnecessary modeling bias. A sufficient condition for unbiasedness is $ p_{\lambda}'(|\theta|)=0 $ for large $|\theta|$. Sparsity: The resulting estimator is a thresholding rule, which automatically sets small estimated coefficient to zero to reduce model complexity. A sufficient condition for sparsity is $\min_{\theta \neq 0}\left\\{|\theta| + p_{\lambda}'(|\theta|)\right\\} > 0 $. Continuity: The resulting estimator is continuous in data $z$ to avoid instability in model prediction. A sufficient and necessary condition for continuity is that the minimum of the function $|\theta| + p_{\lambda}'(|\theta|)$ is attained at 0. * * * It is the last one (continuity) that I don't quite understand how the condition stated results in a continuous estimator. Also Fan mentions: A penalty function satisfying the conditions of sparsity and continuity must be singular at the origin. I also can't grasp why this must be so.


title: Normalizing a Customer Satisfaction Index

text: I'm looking to use the ACSI (American Customer Satisfaction Index) for a customer satisfaction survey. I'm confused by the explanation for how to calculate this. The survey is three questions:   1. What is your overall satisfaction with [our product or service]?   2. To what extent has [our product or service] met your expectations?   3. How well did [our product or service] compare with the ideal [type of offering]?* And the explanation for the calculation is: Some organizations simply normalize and average the three ratings, like this: $$((Satisfaction-1) + (Expectancy-1) + (Performance – 1)) / 27 \cdot 100$$ I don't understand how this works -- why is 1 subtracted from each, and why divide by 27?


title: SARIMA estimation

text: I am trying to manually estimate the non-seasonal components of an SARIMA (p,d,q)x(P,D,Q)[s]. I thought the estimation is going the same way like in ARIMA, but the output says somehow something different. I have an autocorrelation in the acf correlogram and one significance bound at lag 1 in the pacf. That means I have an autocorrelation first order. I'm confused now, why `auto.arima` is giving me the result (0,1,1)x(0,0,1)[12] instead of (1,1,0)x(0,0,1)[12] Here is my code example:               timeseries <- ts(daten, start=c(1955,1), freq=12)          > timeseries           Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec     1955  1.8  1.7  1.5  1.2  1.5  1.5  1.6  1.8  1.5  1.5  1.6  1.3     1956  0.7  0.6  0.4  0.9  0.9  0.8  0.8  0.6  0.6  0.4  0.4  0.2     1957  0.2  0.1  0.6  0.8  0.3  0.4  0.5  0.7  0.8  0.9  1.0  1.3     1958  1.7  1.7  1.4  1.0  0.9  1.3  1.3  1.0  1.5  1.4  1.4  2.2     1959  1.3  1.7  1.7  2.2  2.8  2.5  2.2  2.3  1.8  1.6  1.3  1.4     1960  2.2  1.8  1.9  1.6  1.1  0.8  1.1  1.1  1.1  1.4  1.2  1.2     1961  0.9  1.2  1.3  0.9  0.7  0.8  0.8  1.2  1.0  1.0  1.4  1.0     1962  1.1  0.8  1.1  1.7  2.1  2.0  2.1  2.1  2.0  2.3  2.0  2.3     1963  1.6  1.9  1.6  1.4  1.6  1.8  1.8  1.9  2.5  2.3  2.2  2.1     1964  2.1  2.1  1.9  2.3  2.1  2.0  2.1  1.8  1.0  1.1  1.5  1.4     1965  1.8  1.9  2.0  2.0  2.0  2.0  2.0  2.0  2.7  2.7  3.3  3.1     1966  2.9  3.0  3.3  2.6  3.1  3.4  3.5  3.3  3.0  2.5  1.4  1.1     1967  0.9  1.0  0.4  0.8  0.0  0.0 -0.7 -0.1 -0.5 -0.1  0.3  0.8     1968  0.8  0.5  1.2  1.0  1.2  0.8  1.2  1.0  1.3  1.3  1.6  1.9     1969  2.0  2.2  2.3  2.7  2.4  2.4  2.6  2.5  2.9  2.9  2.8  2.3     1970  2.3  2.5  2.3  2.2  2.2  2.0  1.9  2.2  2.1  2.1  1.9  2.0     1971  1.9  1.8  1.8  1.1  1.6  1.9  1.9   NA           diffts <- diff(timeseries,12)     tsdisplay(diffts, lag.max=36)      ![enter image description here](http://i.stack.imgur.com/2MgzU.jpg) But `auto.arima` is giving me the following output:               auto.arima(timeseries)          Series: timeseries      ARIMA(0,1,1)(0,0,1)[12]                              Coefficients:               ma1     sma1           -0.1280  -0.7260     s.e.   0.0684   0.0584          sigma^2 estimated as 0.07113:  log likelihood=-23.77     AIC=53.54   AICc=53.66   BIC=63.42


> "The set of points in $\mathbb{R}^2$ classified ORANGE corresponds to > {$x:x^Tβ>0.5$}, indicated in Figure 2.1, and the two predicted classes are > separated by the decision boundary {$x:x^Tβ=0.5$}, which is linear in this > case. We see that for these data there are several misclassifications on > both sides of the decision boundary. Perhaps our linear model is too > rigid—or are such errors unavoidable? Remember that these are errors on the > training data itself, and we have not said where the constructed data came > from. Consider the two possible scenarios: > > Scenario 1: The training data in each class were generated from bivariate > Gaussian distributions with uncorrelated components and different means. > > Scenario 2: The training data in each class came from a mixture of 10 low- > variance Gaussian distributions, with individual means themselves > distributed as Gaussian. > > A mixture of Gaussians is best described in terms of the generative model. > One first generates a discrete variable that determines which of the > component Gaussians to use, and then generates an observation from the > chosen density. In the case of one Gaussian per class, we will see in > Chapter 4 that a linear decision boundary is the best one can do, and that > our estimate is almost optimal. The region of overlap is inevitable, and > future data to be predicted will be plagued by this overlap as well. In the > case of mixtures of tightly clustered Gaussians the story is different. A > linear decision boundary is unlikely to be optimal, and in fact is not. The > optimal decision boundary is nonlinear and disjoint, and as such will be > much more difficult to obtain." Can someone please explain to me what the particular scenarios mean? It's from Elements of Statistical Learning by Tibshirani


title: Spearman correlation in the presence of many ties - how to spot a problem?

text: I'm testing the hypothesis that there's a monotonic relationship between two variables. I think I should use a Spearman rank correlation test, since my data don't necessarily meet normality assumptions & have many outliers. However, there are many ties in the independent variable. How can I tell whether the ties are causing me a problem? The data look something like this (R code):               set.seed(0)     x <- rep(1:10, 10)     y <- x + rnorm(length(x), sd=rep(x, 10))      ![enter image description here](http://i.stack.imgur.com/GCwEM.png) One approach I can think of is to add a small random number to each x value many time, and look at the mean/median p-value, like so:               nReps <- 100     pVec <- rep(NA, 100)     for(i in 1:nReps) {       xDodge <- x + rnorm(n=nReps, mean=0, sd=0.0001)       pVec[i] <- cor.test(xDodge, y, method="spearman")$p.value     }     mean(pVec)     sd(pVec)      Does that method seem reasonable? Is there a previously-described method to assess the effect of ties on Spearman's rho, or a similar correlation method that does better with large numbers of ties?


title: Negative density for non-negative variables

text: Having an **integer positive** variable (number of days) in an experiment, I got negative values for the density plots using R. I have read other posts relating to this topic. They admitted that the experimental biases, zero values and values very close to zero are the reasons for this. Can you suggest a link which shows the values close to zero will end up to negative density and why is this so? Next question is that I have only two observations (out of 5000 observations) with zero value. And also there are not many observations with value close to zero as well. However, the density plot in R starts from minus 100. How can I explain it? I will use this link to solve the explanation problem but I also like to know more about it.


title: Test for equality of coefficients from 2 different samples

text: I have a model that I fit on two different samples (each representing a region in a country, with different sample sizes):               region 1: P = a1 + a2*H          region 2: P = b1 + b2*H      which is simply:               region 1: R1 <- lm(P ~ H)          region 2: R2 <- lm(P ~ H)      To test for equality of coefficient estimates, a Wald test would be appropriate. However, the standard Wald test, `waldtest()`, can only be carried out for models applied to the same data set. How can I get around this and test for equality of coefficients when I have 2 different samples? thanks, Diana


title: What does "$\sim$" mean and $A | B \sim C$?

text: I'm not sure if I fully understand the meaning of the symbol, I've seen this symbol in various articles but haven't managed to understand what they implied. I did some reading and it looks like $A \sim B$ means B-Distribution of random variable $A$. Then how would it apply in continuous case as this $$\epsilon_{ij}\sim F$$ $$F\sim DP$$ where $\epsilon$ is noise and DP means Dirichlet process? Then it gets more complicated with $A | B \sim C$ such as Dirichlet process mixture models. Like is the left-hand side denoting $\theta$ given $G$? $$\theta_i|G \sim G$$ $$x_i|\theta_i \sim F(\theta_i)$$ source


title: How to do multiple regression with within-subject independent variables using SPSS?

text: I have two groups of subjects, both undergoing two experimental conditions. During the experimental conditions, the subjects take questionnaires for evaluation of the mood states, and blood pressure is measured at the same time. I want to obtain the correlation between blood pressure and the mood ratings, with simultaneous consideration of the effects of group and condition. Thus, I run a multiple regression, with blood pressure being the dependent variable, and mood ratings, condition as well as group being the independent variables. My question is: Is this a multiple regression with repeated measures? If so, how can I perform the analysis using SPSS. I have tried analysis -> regression -> linear -> set independent variable and dependent variable, but I'm afraid that the settings do not consider a repeated measures design. Thank you very much!!


title: (hierarchical) cluster analysis with non-standard distance

text: My question is triggered by a question that was asked on stackoverflow: http://stackoverflow.com/questions/12198115/using-different-metric-for-hclust- linkage. The thing is this: I can formulate an algorithm for hierarchical clustering with some kind of distance between original objects, and I'll also need a distance between objects and clusters of objects (linkage). Let's say the linkage attaches a "handle" to the cluster which is a (hypothetic) object that can be used with the distance function to calculate the distance from another object (or cluster, if it has such a "handle"). For consistency, I'd use _the same_ distance function for the initial calculation of distance between single objects and later for distance between objects or "cluster handles". So far, so good. What I'm wondering about is: what happens if the cluster algorithm is fed with a distance matrix instead of the data matrix (as it is e.g. the case with R's `hclust`).   * How do I know which distance function on the _distance matrix_ agrees with what distance function on the _data matrix_? Examples?    * The more basic maths/statistics question underlying this is: please explain to me (chemist, no mathematician/statistician) what the meaning of e.g. Euclidean distance between rows of a distance matrix is?


title: How to select significant variables for multi-factor regression?

text: I have a set of data $Y \in \mathbb{R}^n$ and a set of $k$ factors $\\{F_1, ..., F_K\\}$ with $F_i \in \mathbb{R}^n ~ \forall i$. I would like to perform a factor analysis which consits in finding $\beta \in \mathbb{R}^k$ such that: $$Y_t = \alpha + \sum_{i=1}^k \beta_i {F_i}_t + \varepsilon_t$$ I believe that I would need to make sure that the factors are independent from each other, and I wonder what process (statistical test) I should use to determine which factors I shall remove from the framework because they are too correlated. There are two main questions:   1. How do I find and remove the factors that are not linearly independent?   2. How do I determine which of the remaining factors are significant?


title: How to compare parameters between different types of generalized linear model?

text: I am doing some work on the effects of collinearity on different types of model (OLS, binomial logistic, ordinal logistic, multinomial logistic and maybe others). I have found the perturb package in R, which will be very useful for this. This package shows how small changes in the data affect the parameter estimates of a model. The perturb package perturbs the data that is used to build the model, and then re-runs the model - a little like a bootstrap. In this particular case, I chose to perturb the independent variables, but I could also perturb the dependent variable. However, I am struggling to figure out the best way to compare these effects. For example, I ran perturb on an OLS model and got               ----------                Impact of perturbations on coefficients:                           mean  s.d.    min   max              x0          0.014 0.264 -0.743 0.604              x1          1.451 0.347  0.771 3.772              x2          0.847 0.321  0.308 2.817              x3          2.222 0.243  1.129 2.746         ----------      Then I dichotomized the DV and ran perturb on a logistic model and got               ----------                Impact of perturbations on coefficients:                           mean  s.d.    min    max              x0          0.119 0.234 -0.477  0.732              x1          0.918 0.819 -6.612  2.215              x2          0.247 1.008 -9.413  1.175              x3          1.533 2.292  0.574 24.083     ----------      The latter certainly _looks_ like it was affected much more, but I don't know of any more formal ways to compare these parameter estimates and, especially, the size of the sd of the different estimates. Here is a link to the perturb package and to the R version Thanks


title: What is the advantage of having balanced panel data rather than unbalanced?

text: I studied the standard econometrics textbooks about panel data, but most textbooks only mention the difference between balanced and unbalanced panels. The advantage of having balanced panel data is not usually explained. I would like to know: **What is the advantage of having a balanced panel?** I believe unbalanced panels are much more common in real research.


title: Size of the data set for measuring variation

text: This might be a simple question but I am unable to find a solution on the web. I would like to measure variation of some numerical data. What is the minimum number of values I should use for measuring their standard deviation or coefficient of variation? Basically my plan is to derive groups from a set of performance metric values (all are positive integers). I would like to choose "variation" measure for grouping. My idea is as follows: For example, I have 16 positive integer numbers. First I would like to sort them then divide them as four groups. Each group contains 4 integers and I will calculate variation in terms of either standard deviation or coefficient of variation. If the difference in the variation of two groups is less than a threshold then I will merge them, and etc. Here my question is, can I take a data size of 4 for calculating variation? Please help me.


title: How to analyse multiple between subjects factors in a mixed ANOVA in SPSS?

text: I am performing a Mixed between subjects Anova where I have multiple time points that subjects have completed a measure at.   * The within subjects factor is time.    * The between subjects factors include gender, ethnicity, refugee status and English language.  My two questions are   1. **Combined or separate ANOVAs** : Is it possible in SPSS to do this as one analysis or do I need to do multiple Anovas one for each of my between subjects factors? If it is possible to do it in one model can you explain how specifically in SPSS?   2. **Fixed versus random factors** : My advisor told me to be sure to know whether I am using the models as fixed effects or random effects but she herself doesnt really understand this. I understand the difference ie fixed effects are like gender versus random effects that are random sample of the total # of levels in the population however I am not sure how to deal with this in SPSS or even tell what the default is. I performed my analyses using Analyze, General Linear Model, Repeated Measures.


title: Observation (evidence) in dynamic Bayesian networks

text: I am wondering how the probabilities of the observation nodes in dynamic Bayesian networks are set. I want to know whether the probabilities are monitored or are given by sensors? So, what does the term `observation nodes` actually mean, and from what can I compute their probabilities?


title: What processes could generate Laplace-distributed (double exponential) data or parameters?

text: Lots of distributions have "origin myths", or examples of physical processes that they describe well:   * You can get normally distributed data from sums of uncorrelated errors via the Central Limit Theorem   * You can get binomially distributed data from independent coin flips, or Poisson-distributed variables from a limit of that process   * You can get exponentially distributed data from waiting times under a constant decay rate. And so on. But what about the Laplace distribution? It's useful for L1 regularization and LAD regression, but it's hard for me to think of a situation where one should actually expect to see it in nature. Diffusion would be Gaussian, and all the examples I can think of with exponential distributions (e.g. waiting times) involve non-negative values. Thanks in advance for your thoughts.


title: Re-parametrizing Gaussian distribution

text: Let there be observed data points $X = {X_1.. X_n .. X_N}$, where each $X_n \in R^D$. Lets assume these are distributed as a Gaussian $X \cong \mathcal N(\mu,\Sigma)$. Let us also assume that the mean has a Normal prior and covariance has Inverse Wishart prior. Since these are conjugate priors, given observed data $X$, we can draw posterior samples for mean and covariance using standard posterior update equations. What we have so far is just the usual Bayesian stuff. Let the data points be formed from a union of two sets of data $X = {X^+} \cup {X^-}$. The data is still distributed as Gaussian. However lets assume that while drawing posteriors for mean and covariance, we have to give more weight to data points in ${X^+}$ than ${X^-}$, with the weight itself being some (hyper)parameter. Note that we are _NOT_ talking about mixture of gaussians. We dont say there are 2 mixtures and the data might have arised from one of these mixtures. We say the data still arises from ONE gaussian, with the mean and covariance primarily according to ${X^+}$, but perhaps altered(shifted?) by some proportion because of ${X^-}$ influence. I am not sure how i can write the distribution for such a thing. I initially thought i can write $X \cong \mathcal N(A\mu^+ + B\mu^-,A\Sigma^+ + B\Sigma^-)$. We then draw $\mu^+, \Sigma^+$ from ${X^+}$ and similarly $\mu^-, \Sigma^-$ from ${X^-}$ by usual methods. We treat $A$ and $B$ as some (hyper)parameters. If i write like this, what is my assumption of the data. Am i implying $X = {X^+} + {X^-}$ rather than a union? Or the implication is too non-sensical to even worth a description? Why would i want to do this - thats a lengthy explanation and would distract the specific question. The specific question is - "I have data partitioned into two sets. The data arises from a gaussian distribution. I want to draw posterior samples for the gaussian parameters according to one of the set. However the other set also has (marginal) influence on the parameters. How can I re-parametrize the distribution to account for this?" Ps: Any approximate re-parametrization is fine. I've used the word "sampling", hopefully that wouldn't distract the specific question to ideas such as lets reject samples based on some function. Any references to the literature that does something similar is also welcome. Thanks in advance for suggestions.


title: How many factors to retain in EFA (Exploratory Factor Analysis)?

text: I have multiple measures from a scale and I want to determine the best factorial structure using EFA, in SPSS. I don't understand if I can retain the number of factors that I want to retain because they seem to be the most theoretically valid, or if I must keep the number of factors the program gives me based on Kaiser criterion or a scree plot. Say, I want to extract a total of 4 factors, after eliminating cross-loadings and items with small communality scores. Based on these criteria, the program extracts 6 (hypothetical) factors. Can I say I used EFA and retained 4 factors? Is that correct? If I'm not clear, please state it. Any help greatly appreciated!


title: Testing difference of means between two samples, without distribution

text: I have a data set where I am trying to compare a difference of means between two different samples. However, I do not have the standard deviation of the mean for each sample. The sample data consists of a date, a number of visits, and a total number of times the event of interest. For each visit, the event of interest can be triggered multiple times. Example data below.               +----------+-----------+-------+--------+------------+     |   Date   | Treatment | Vists | Event1 | Avg.Event1 |     +----------+-----------+-------+--------+------------+     | 1/1/2012 | A         |  5392 |  12390 | 2.297      |     | 1/1/2012 | B         |  5489 |  11499 | 2.094      |     | 1/2/2012 | A         |  6030 |  12332 | 2.045      |     | 1/2/2012 | B         |  6200 |  12003 | 1.935      |     +----------+-----------+-------+--------+------------+      So with a data set like that, I assume I can't use hypothesis tests directly to compare the sample means. Is there a statistical method to test the difference between treatment A and treatment B? Any suggestions or pointers in the right direction would be greatly appreciated. Thank you.


title: Interaction variables work differently when population is split

text: Let's say I have two predictors to predict financial risk: Gender and shopping habits. Gender has levels of "Male" and "Female", while shopping habits has "quick shopper" and "slow shopper". I am quite surprised to see this: When I split the population into "male" and "female", shopping habits as a single variable work quite well to separate the risk on female population (not for males). This also intuitively makes sense. But when I create an interaction variable of "gender:shopping habits" and look at logistic regression, there is no risk separation whatsoever even on the dummy created "female:quick shopper". Could there be an explanation of this? Note: I should add that the initial working performance is cross-validated (not an overfitting result). Note2: I edited the question for better understanding.


title: time series (asynchronous/variable) sample rate conversion

text: I am having a problen in preparing my data. The data is a GPS track recorded from a lengthy car-trip. This is recorded with a constant frequency - I have a datapoint every 200ms. Now I have calculated the distance that was travelled up to each point. What I want now, is to interprete the distance as time base of my series. This is of course strictly increasing, but not with an constant frequncy. How far it jumps with each sample is dependend on the velocity. So I want to end up with a series of data-point that haven an equal distance from each other. Other attributes, like speed or heading, needs to be mapped to that form, which implies to sometimes insert new data points and sometimes remove some. Now is there an algorithm out there that can accomplish this? I need to implement this as a KNIME-node, but I have access to R, Python, Java, (Matlab). Hints appreciated, daniel


title: Who said, "let the data speak for themselves"?

text: Who said: > "Let the Data Speak for Themselves" \-- Ronald Fisher or John Tukey?


title: K-Means clustering after first iteration

text: In _k-means clustering_ we initially pick $k$ random centroids and assign the given data to one of these $k$ centroids (which ever is nearest). After this we create new centroids by taking the mean of the assigned points. However there might be case that the initially selected random centroids may not be nearest to any point in the dataset and hence no points would be assigned to these centroids. So in such case what should be done in the step of creating new centroids?


title: Simple condition for asymptotic normality of MLE

text: Are these conditions sufficient for asymptotic normality of MLE? If it is, pls. let me know the references. 1-First and second derivatives of $\ell(\theta,\eta)$ are defined. 2- The Fisher information matrix be non-singular and continuous with respect to the parameters $\theta$ and $\eta$.


title: Simple linear regression between continuous outcome and dichotomous predictor

text: Can I run simple linear regression on SPSS between continuous outcome such as exam results and a dichotomous predictor such as smoking or not (0 and 1).


title: Binomial random variable, solving for number of trials.

text: Assuming that each flight is independent of the others, and each flight has a probability of $10^{-9}$ of crashing independent of other flights, how many flights are needed for the probability of the aircraft crashing at least once to reach $0.01$. I'm trying to use this: $0.01=\left(\array{n\\\1}\right)(10^{-9})(1-10^{-9})^{(n-1)}$ when I solve for n using the definition of the binomial coefficient I'm getting a negative value.


title: Temporal autocorrelation in perMANOVA?

text: I have a data set where samples are collected once per year for 15 years at a number of sites. I am worried that these data are temporally autocorrelated and was trying to figure out if I need to address that. However, the only time I will be using degrees of freedom with these data is in a perMANOVA. This test calculates a pseudo F-statistic by permuting the rows. I can't figure out if the exchangebility assumption means that I don't need to worry about autocorrelation at all (i.e., permuting rows will simply destroy the temporal structure, which I am not interested in anyway) or if it means that I can't use a perMANOVA even if I accounted for autocorrelation? Edit: I am editing this in the hopes that clarification will help get it answered. The perMANOVA user's guide says: _"Recall that for traditional one-way ANOVA, the assumptions are that the errors are independent, that they are normally distributed with a mean of zero and a common variance, and that the treatment effects are additive. In the case of a one-way analysis, the PERMANOVA test using permutations assumes only that the observation units are exchangeable under a true null hypothesis. There are no explicit assumptions regarding the distributions of the original variables; they are certainly not assumed to be normally distributed. However, implicit in the notion of exchangeability is the notion of independence, for if observations are correlated with one another (e.g., temporally or spatially), then randomly shuffling them will destroy this kind of inherent structure, if it is there. Thus, in general, we would assume that the observation units are independent of one another."_ The meaning of this is ambiguous to me for the reasons stated in the first paragraph. I can't find any techniques for testing/correcting autocorrelation with perMANOVA, which maybe means that it isn't a problem to worry about? User's guide: https://www.stat.auckland.ac.nz/~mja/prog/PERMANOVA_UserNotes.pdf Original 2001 paper describing technique: http://50-56-175-137.static.rackspacecloud.com/PDF/MUVE/6_NewMethod_MANOVA1_2.pdf


title: Rescaling input features for neural networks regression

text: In Neural Nets for the regression problem, we rescale the continuous labels consistently with the output activation function, i.e. normalize them if the logistic sigmoid is used, or adjusted normalize them if tanh is used. At the end we can restore original range but renormalizing the output neurons back. Should we also normalize input features? And how? For example, if hidden activation differs from the output activation? E.g. if hidden activation is TANH and output activation is LOGISTIC, should the input features be normalized to lie in [0,1] or [-1,1] interval?


title: Testing on non-normal distributed discrete values

text: I have the following observations Oberservation ; Count -1.67 ; 726  18.33 ; 33 148.33 ; 15 This is obviusly not normal distributed :S How can I make a test for $H_0: \mu = 0$ or even better is it possible to make a confidence interval for the mean?


title: polynomial regression model

text: I have 8 parameters with 3 levels and the number of tests required to study this case is 27. Can I generate polynomial regression model to describe interaction of these parameters?


title: Method of Composition to sample from a t density

text: I got stuck with this, I will appreciate a lot any help. I need to make an R program in order to run this algorithm (in the photo below), with simulated data. The question is to use the method of composition to sample from a t-distribution via sampling normal distributions. ![this is the algorithm](http://i.stack.imgur.com/rksyP.png) I have been working with pre-existing code from a textbook in Bayesian Statistics, as follows:               # SAMPLING FROM NONINFORMATIVE PRIOR     attach(ALP)     objects(2)          alp <- alkfos[artikel==0]     talp <- 100*alp^(-1/2)          mtalp <- mean(talp)     vartalp <- var(talp)     sdtalp <- sqrt(var(talp))     ntalp <- length(talp)          # sample from sigma then from mu given sigma  (1000 times),      #  make use of theoretical results     # Figure 4.3     set.seed(777)          unl.mu <- 6.8     upl.mu <- 7.4     unl.sigma2 <- 1.2     upl.sigma2 <- 2.5          nmu <- 50     nsigma2 <- 50     ntilde <- 50      .               mu <- seq(unl.mu,upl.mu,length=nmu)     mus <- (mu-mtalp)/(sdtalp/sqrt(ntalp))     pmu <- dt(mus,ntalp-1)*sqrt(ntalp)/sdtalp     sigma2 <- seq(unl.sigma2,upl.sigma2,length=nsigma2)     sigma2s <- (ntalp-1)*vartalp/sigma2     psigma2 <- dchisq(sigma2s,ntalp-1)*(ntalp-1)*vartalp/(sigma2^2)          par(mfrow=c(2,2))     par(pty="s")     par(mar=c(5,4,4,2)+0.1)          nsimul <- 1000      .               #posterior distribution of sigma^2         phi.u <- rchisq(nsimul,ntalp-1)     theta.u <- 1/phi.u     sigma.u <- sqrt((ntalp-1)*vartalp*theta.u)          # posterior summary statistics     meansamplesigma2 <- mean(sigma.u^2)     sdsamplesigma2 <- sqrt(var(sigma.u^2))     unl.conf <- meansamplesigma2 -1.96*sdsamplesigma2/sqrt(nsimul)     upl.conf <- meansamplesigma2 +1.96*sdsamplesigma2/sqrt(nsimul)     meansamplesigma2;sdsamplesigma2;unl.conf;upl.conf      .               quantile(sigma.u^2, c(.025,.975))            hist(sigma.u^2,xlab=expression(sigma^{2}),     ylab="", main="", nclas=30,prob=T,cex.lab=1.5,col="light blue",          cex.axis=1.3)     lines(sigma2,psigma2,lwd=3,lty=1,col="blue")     text(1.5,2.3,"(a)",cex=1.2)          #conditional posterior distribution of mu given sigma^2     mucond1.u <- rnorm(nsimul,mtalp,sdtalp/sqrt(ntalp))     mucond2.u <- rnorm(nsimul,mtalp,2.5/sqrt(ntalp))          #hist(mucond1.u,xlab="MU values, given SIGMA = SAMPLE VARIANCE",     #ylab="HISTOGRAM",main="(b)",nclas=30,prob=T)     #hist(mucond2.u,xlab="Posterior sampled MU values,given SIGMA = 2.5",nclas=30,prob=T)      .               #marginal posterior distribution of mu given sigma^2     mu.u <- rnorm(nsimul,mtalp,sigma.u/sqrt(ntalp))          # posterior summary statistics     meansamplemu <- mean(mu.u)     sdsamplemu <- sqrt(var(mu.u))     unl.conf <- meansamplemu -1.96*sdsamplemu/sqrt(nsimul)     upl.conf <- meansamplemu +1.96*sdsamplemu/sqrt(nsimul)     meansamplemu;sdsamplemu;unl.conf;upl.conf          quantile(mu.u, c(.025,.975))            hist(mu.u,xlab=expression(paste(mu)),ylab="",main="",     nclas=30,prob=T,cex.lab=1.5,col="light blue",          cex.axis=1.3)     lines(mu,pmu,lwd=3,lty=1,col="blue")     text(6.9,5,"(b)",cex=1.2)          #joint posterior distribution of mu and sigma^2     plot(mu.u,sigma.u^2,xlab=expression(paste(mu)),ylab=expression(paste(sigma^{2})),          main="",cex.lab=1.5,col="blue",          cex.axis=1.3)     text(6.9,2.4,"(c)",cex=1.2)      .               #sample from ppd          future.u <- rnorm(nsimul,mu.u,sigma.u)          hist(future.u,xlab=expression(tilde(y)),nclas=30,prob=T,main="",          cex.lab=1.5,ylab="",ylim=c(0,0.3),col="light blue",          cex.axis=1.3)     text(3.5,0.3,"(d)",cex=1.2)          py.unl <- mtalp -qt(0.999,ntalp-1)*sdtalp*sqrt(1+1/ntalp)     py.upl <- mtalp +qt(0.999,ntalp-1)*sdtalp*sqrt(1+1/ntalp)          tildey <- seq(py.unl,py.upl,length=ntilde)     tildeys <- (tildey-mtalp)/(sdtalp/sqrt(1+1/ntalp))     ptildey <- dt(tildeys,ntalp-1)*sqrt(1+1/ntalp)/sdtalp          lines(tildey,ptildey,lwd=3,lty=1,col="blue")      I just have changed the first part of the code generating a t-distribution as follows:               x<-rt(10000, 1)     mx<-mean(x)     varx<-var(x)     sdx<-sqrt(var(x))     nx<-length(x)      And then, changing the name of the variables in the corresponding parts (also erasing those lines where the original code makes some transformations). Doing that, I get a result, but, I am not sure if that is the way. Thanks!


title: Cluster analysis in R produces reversals on dendrogram

text: I'm attempting to perform hierarchical agglomerative cluster analysis in R. However, when I use particular clustering methods, I get reversals (upward branching) in the resulting tree, which violates the ultrametric property. ![enter image description here](http://i.stack.imgur.com/WSPo2.jpg) The two methods are: UPGMC and WPGMC (methods="median" and "centroid" in `hclust`). Legendre & Legendre in their Numerical Ecology book suggest some reasons why this may occur (Section 8.6). However, they provide no solutions to rectify the issue and convert the trees to ultrametric. I'm curious: is this an unavoidable consequence of the data and the clustering method, or is there a way that I can produce a tree that satisfies the ultrametric property using these two methods? Here is an example data set and R code to play with:               #Generate data frame with mixed continuous and categorical trait data for 10 species     set.seed(91)     (df=data.frame(trait1=runif(10,0,10),trait2=runif(10,0,10),                    trait3=sample(letters[1:3],10,replace=T),row.names=paste("sp",1:10,sep="")))          #Generate Gower dissimilarity matrix from trait data     library(cluster)     (dist.gower=daisy(df,metric="gower"))          #Create a vector of clustering methods     tree.methods=c("ward","single","complete","average","mcquitty","median","centroid")       #Build the trees using each method     trees=lapply(tree.methods,function(i) hclust(dist.gower,method=i))       #Plot the trees     par(mfrow=c(4,2))     for(i in 1:length(trees)) {plot(trees[[i]])}     #The last two trees have reversals...cannot be converted to ultrametric!


title: GLM post hoc with non-parametric tests?

text: I have a question regarding the appropriate use of comparisons for independent samples (3 factor levels). Overall sample size is N=546, subsamples: 218 or 228 or 100), convenience sampling, stratified. I use ANOVA with post hoc Tukey, if Levene indicates variance homogeneity. **My question: If there no equal variance: What to do?** a. Using Kruskal-Wallis? For post hoc: Do I use Mann-Whitney for pairwise comparisons (if so: what about the alpha inflation?)? b. Still using GLM but with non-parametric post hocs? If so: which one? Tamhane? Dunnett or Games-Howell? What about F, can I still use it? Isn’t H then the more appropriate statistic? Why gives GLM anyway the option of running on-parametric tests of homoscedasticity is an assumption (I am confused) Thanks for any input on this one


title: $2^{5-2}$ Design : Alias Structure

text: $2^{5-2}$ Design Design Generators: $D=AB\quad E=AC$ Defining Relation: $I=ABD=ACE=BCDE$ $$\text{Aliases}$$ $$A=BD=CE=ABDE$$ $$B=AD=ABCE=CDE$$ $$C=ABCD=AE=BDE$$ $$D=AB=ACDE=BCE$$ $$E=ABDE=AC=BCD$$ $$BC=ACD=ABE=DE$$ $$CD=ABC=ADE=BE$$ But the alias structure in the book `Design and Analysis of Experiments` by Douglas C. Montgomery [Appendix XII (c)] is as following: $$\text{Aliases}$$ $$A=BD=CE$$ $$B=AD=CDE$$ $$C=AE=BDE$$ $$D=AB=BCE$$ $$E=AC=BCD$$ $$BC=DE=ACD=ABE$$ $$CD=BE=ABC=ADE$$ Instead of $A+BD+CE+ABDE$, this is $A+BD+CE$ in the book and so on.   * Why ?


title: Predicting with both continuous and categorical features

text: Some predictive modeling techniques are more designed for handling continuous predictors, while others are better for handling categorical or discrete variables. Of course there exist techniques to transform one type to another (discretization, dummy variables, etc.). However, are there any predictive modeling techniques that designed to handle both types of input at the same time without simply transforming the type of the features? If so, do these modeling techniques tend to work better on data for which they are a more natural fit? The closest thing that I know of would be that usually decision trees handle discrete data well and they handle continuous data without requiring an _up front_ discretization. However, this isn't quite what I was looking for since effectively the splits on continuous features are just a form of dynamic discretization. For reference, here are some related, non-duplicate questions:   * How should decision tree splits be implemented when predicting continuous variables?   * Can I use multiple regression when I have mixed categorical and continuous predictors?   * Does it ever make sense to treat categorical data as continuous?   * Continuous and Categorical variable data analysis


title: Difference in Effect for Time-Varying Covariate in Survival Analysis

text: I am analyzing offender recidivism data using survival analysis. In particular, I am looking at the risk of getting arrested as a function of employment. When I treat employment as a time-varying covariate, it is associated with a 21% reduction in the risk of recidivating. However, if employment is included without reference to time {i.e., 0/1 indicator of getting a job during the time at risk that does not change over time} employment is associated with a 75% reduction in the hazard of recidivism. Can someone explain why I'm seeing this difference?


title: Weights of random sets of random 32-bit strings

text: I have random sets of $N$ random 32-bit strings, where all bits are i.i.d. with $\mathbb{P}(0) = \mathbb{P}(1) = 1/2$. Define   $\ \ \ \ $weight( 32-bit x ) = number of 1 bits in x, i.e. Hamming distance to 0   $\ \ \ \ $minweight( set $X$ ) = min$_{\text{x} \in X}$ weight( x ) How does the average minweight( $X$ ) vary with $N$, the size of $X$ ? My $N$ is small enough that it doesn't matter whether $X$ has duplicates, or not.   A back-of-the envelope estimate would be fine. (A harder question is to analyze the algorithm for bit-string-nearest- neighbour-searching on stackoverflow, but that seems rather specialized. If anyone could suggest another forum or reference, I'd appreciate it.)


title: ANCOVA in R and setting contrasts

text: Using the following data, I need to determine if, independent of mass, the observed response is associated with sex, location, and/or age               example<- data.frame(response=c(2.401,2.588,2.293,2.880,2.655,2.830,3.165,2.665,3.126,2.973,1.725,1.889,2.631,1.750,2.271,2.347,2.173,1.962,1.599,2.297,1.894,2.030),          mass=c(460,510,500,475,450,470,600,570,480,560,410,370,425,385,400,420,490,450,390,465,360,490),      loc=c("K","K","K","K","K","K","K","K","K","K","R","R","R","R","R","R","R","R","R","R","R","R"),     sex=c("F","M","M","F","M","F","M","M","F","M","F","F","M","F","F","F","M","F","F","M","F","F"),     age = c(48,49,47,45,49,49,50,48,49,48,47,48,47,50,48,50,48,50,48,50,50,50))      When running an ANCOVA in R, I've noticed that some examples set orthagonal contrasts and other do not. Do I need to set the contrasts, and how would I code it for this data set? I'm guessing everything else is straightforward (with pkg 'car' installed):               model.1=aov(response~mass+loc+sex+age, data=example)     Anova(model.1, type="III")      Summary.lm(model.1)


title: Collaborative Filtering, matrix decomposition, and incorporate various kinds of data

text: Just about every matrix factorization (e.g., SVD++) has some matrix that includes a n users and m (e.g.,) ratings. Here is my question, how do you include information like demographic information, information on what the user has said about a product (like, "it's cheesy", or "awesome!"). What I'm really failing to see is how to incorporate this extra info into a matrix factorization approach.


title: How to do logistic regression in R when outcome is fractional?

text: I'm reviewing a paper which has the following biological experiment. A device is used to expose cells to varying amounts of fluid shear stress. As greater shear stress is applied to the cells, more of them start to detach from the substrate. At each level of shear stress, they count the cells that remain attached, and since they know the total number of cells that were attached at the beginning, they can calculate a fractional attachment (or detachment). If you plot the adherent fraction vs. shear stress, the result is a logistic curve. In theory, each individual cell is a single observation, but obviously there are thousands or tens of thousand of cells, so the data set would be gigantic, if it was set up in the usual way (with each row being an observation). So, naturally, my question (as stated in the title) should make sense now. How do we do a logistic regression using the fractional outcome as the D.V.? Is there some automatic transform that can be done in glm? Along the same lines, if there were potentially 3 or more (fractional) measurements, how would one do this for a multinomial logistic regression?


title: Modify Levenberg Marquardt

text: Could you kindly let me know if it is possible to print parameters to the console as they are being optimized (using the levenberg marquardt algorithm) in python/matlab? I have been trying to do so but to no avail. I would need to modify the printed parameters before supplying them into the algorithm before proceeding with the optimization. Again, I also wish to know if anyone out here has any idea on how to set a set of "seeded" initial values before starting the Levenberg Marquardt algorithm in either python/matlab? This is vital, given the fact that the choice of initial parameters play a key role in the optimization procedure. Thanks alot.


title: Gaussian Mixture Model parameters from density

text: How do I estimate parameters of subpopulations in a 1D gaussian mixture model when I already have density (measured on a grid) of the mixture? All the algorithms I can find (like the well-known EM algorithm) assume that only samples from the mixture are available. My experiment directly yields density values, not separate samples. Of course I could sample manually from that density, but I feel there must be a better way


title: What can I say about my logistic regression output?

text: I'm analyzing customer satisfaction for a certain program. My response variable is ordinal- where 0 means dissatisfied and 1 means satisfied. I have 10 predictor variables. I used Binary logistic Regression on my data, and the results showed that 6 of my 10 variables are statistically significant. My question is, what can I safely say about my output? Can I say that there is significant evidence to suggest the variables have an affect on the response?


title: Fitting a correlation function to my data

text: I have some spatial data and I want to fit a correlation function may be exponential or gaussian to my data. I could calculate the correlation for different pairs of my spatial point to get some empirical correlation function. What is the best way to fit some standard correlation functions to my empirical one. I tried to use the gstat package to estimate the variogram. However, I had some issues. I generated some synthetic data with the given correlation matrix. I mean I generated some samples from a multivariate gaussian distribution with the given correlation matrix. I am assuming it to be homogeneous. My first row of correlation matrix is something like               1   0.683429646106259   0.437789623184770   0.262915008717623   0.196106083531309   0.149045509783630   0.150509368156993   0.0707581842359845  0.0946249108754618  0.0778169918490679      So it kind of decreases with distance. The range is somewhere around at distance 10 when the correlation reaches approx 0.05. Now when I tried to use the gstat package. I had to use the following data.               "x" "y" "data"     1 1 0.44168     1 2 0.099929     1 3 0.026909     1 4 -0.041947     1 5 -0.24288     1 6 0.78454     1 7 0.16761     1 8 -0.086277     1 9 0.023661     1 10 -0.26036     2 1 0.71725     2 2 -0.14143     2 3 -1.0092     2 4 -0.28496     2 5 -0.17743     2 6 0.22874     2 7 -0.94645     2 8 -0.14165     2 9 0.23027     2 10 0.50463     3 1 -1.1799     3 2 0.31035     3 3 0.83902     3 4 0.057151     3 5 -0.77414     3 6 -0.57698     3 7 0.034864     3 8 0.37105     3 9 -0.017573     3 10 -0.80297     4 1 -1.1183     4 2 0.04495     4 3 -0.43193     4 4 -0.32199     4 5 -0.30001     4 6 0.032265     4 7 0.76654     4 8 0.047539     4 9 -0.70933     4 10 -0.39096     5 1 0.8237     5 2 0.48326     5 3 -0.18275     5 4 -0.2011     5 5 -0.45547     5 6 0.056931     5 7 0.38949     5 8 -0.29351     5 9 -0.82567     5 10 -0.6183     6 1 1.2794     6 2 0.55406     6 3 1.0454     6 4 -0.80348     6 5 0.10936     6 6 -1.5121     6 7 0.86288     6 8 0.18994     6 9 -0.019555     6 10 -0.68697     7 1 1.7539     7 2 1.2916     7 3 0.41383     7 4 -0.24399     7 5 -0.023723     7 6 -1.0232     7 7 0.92783     7 8 0.64965     7 9 0.95411     7 10 0.19481     8 1 0.67385     8 2 1.1123     8 3 0.96762     8 4 0.69418     8 5 -0.025346     8 6 -0.69347     8 7 0.78425     8 8 1.082     8 9 1.1929     8 10 1.3393     9 1 0.87725     9 2 0.30799     9 3 -0.14401     9 4 -0.24109     9 5 0.3923     9 6 1.1385     9 7 1.1368     9 8 0.13259     9 9 0.86029     9 10 0.62073     10 1 -0.043159     10 2 0.13516     10 3 0.44809     10 4 0.1643     10 5 0.66301     10 6 1.9096     10 7 1.1722     10 8 1.1109     10 9 1.1154     10 10 1.0472      This is from just one sample. I thought of providing it more samples but I had to repeat the coordinates and it gave me rubbish variogram like this ![enter image description \[here](http://ifasstat.ifas.ufl.edu/STA%206934SP10/notes/Fitting%20a%20model%20variogram.pdf) However, I tried with just one sample and this is what I got for fitting a model variogram                 model      psill    range     1   Nug 0.06875421 0.000000     2   Exp 0.48353408 1.931677      ![enter image description here](http://i.stack.imgur.com/e8As0.png) ![enter image description here](http://i.stack.imgur.com/NpwOz.png) I definitely cannot take these params for correlogram. In the case of correlation function, range is the distance at which the correlation is 0.05. If I look at my actual correlation function its approx at distance 10. So it doesn't match at all. Is it something wrong with the way I calculated the variogram or something else. I followed the tutorial here I am not sure but I think I am not supposed to get the nugget either. I added no noise when I sampled from the gaussian. Is it because I have less samples? When I try to use all the samples, I get something like this as my empirical variogram which is not correct. Here is my data               "x" "y" "aod"     1 1 -0.75488     1 2 -0.78846     1 3 -0.96701     1 4 0.21699     1 5 0.52046     1 6 0.28913     1 7 0.12993     1 8 0.37072     1 9 0.58311     1 10 -0.16038     2 1 -0.2117     2 2 -1.2046     2 3 -1.0581     2 4 -0.67419     2 5 1.0385     2 6 1.1686     2 7 0.14002     2 8 -0.16956     2 9 0.25615     2 10 0.21616     3 1 1.0365     3 2 -1.3079     3 3 -0.51374     3 4 -0.35841     3 5 0.25237     3 6 0.97418     3 7 -0.10198     3 8 -0.21238     3 9 -0.39775     3 10 -0.93731     4 1 0.50214     4 2 -0.21804     4 3 -0.63425     4 4 0.68123     4 5 1.695     4 6 1.6203     4 7 -0.0013243     4 8 0.15842     4 9 0.29372     4 10 -1.7714     5 1 0.66179     5 2 -1.0807     5 3 -0.65248     5 4 0.15395     5 5 0.30281     5 6 -0.075867     5 7 0.31152     5 8 -1.2481     5 9 -0.8156     5 10 -1.3282     6 1 0.20887     6 2 -0.316     6 3 0.62447     6 4 -0.073839     6 5 -0.34457     6 6 -0.047924     6 7 0.044455     6 8 0.40695     6 9 -0.030344     6 10 -0.078803     7 1 -0.43136     7 2 0.46374     7 3 0.76833     7 4 1.1933     7 5 -0.72779     7 6 0.68372     7 7 0.3142     7 8 -0.78668     7 9 -0.41529     7 10 -1.1691     8 1 0.12466     8 2 0.26433     8 3 0.68538     8 4 0.55875     8 5 0.52583     8 6 -0.127     8 7 0.050293     8 8 -0.18494     8 9 0.25105     8 10 -0.84661     9 1 -0.22072     9 2 -0.11841     9 3 1.0484     9 4 1.2102     9 5 0.77315     9 6 0.76342     9 7 0.45916     9 8 0.0039676     9 9 0.014329     9 10 -0.64462     10 1 0.80869     10 2 -0.35877     10 3 0.73614     10 4 -0.4413     10 5 0.73746     10 6 0.31207     10 7 0.14122     10 8 1.3785     10 9 -0.89216     10 10 0.095206     1 1 1.2769     1 2 1.0405     1 3 -0.0031401     1 4 0.047758     1 5 -0.74042     1 6 -0.69732     1 7 -0.041547     1 8 -0.22316     1 9 0.066819     1 10 -0.010506     2 1 0.96758     2 2 0.31139     2 3 -0.17409     2 4 0.25964     2 5 -0.85022     2 6 -0.71892     2 7 -0.75265     2 8 -1.1433     2 9 -1.4629     2 10 -1.5973     3 1 0.092515     3 2 -0.22445     3 3 0.15134     3 4 -0.46512     3 5 -0.06476     3 6 -0.23128     3 7 -1.3388     3 8 -0.97593     3 9 -1.448     3 10 -0.43345     4 1 0.61175     4 2 0.61255     4 3 -0.10015     4 4 0.33106     4 5 0.34014     4 6 -0.23341     4 7 -0.86573     4 8 -1.2642     4 9 -0.62796     4 10 -1.1511     5 1 0.31983     5 2 0.53389     5 3 -0.76606     5 4 -0.73465     5 5 0.84788     5 6 0.47175     5 7 0.53173     5 8 -0.11972     5 9 0.21588     5 10 -0.46452     6 1 -0.7283     6 2 -0.70479     6 3 -1.1243     6 4 -0.15782     6 5 0.17425     6 6 -0.26241     6 7 1.1181     6 8 -0.15918     6 9 -0.69752     6 10 -0.68888     7 1 -0.56861     7 2 -1.2355     7 3 -0.6948     7 4 0.58936     7 5 0.056865     7 6 0.64748     7 7 -0.58011     7 8 0.070476     7 9 -0.92674     7 10 0.24241     8 1 -0.46411     8 2 -0.44081     8 3 0.54728     8 4 1.5323     8 5 -0.05007     8 6 0.98995     8 7 -0.183     8 8 -0.91211     8 9 0.00081732     8 10 -0.18632     9 1 -0.032806     9 2 0.52061     9 3 1.4889     9 4 0.56336     9 5 -0.18504     9 6 -0.79679     9 7 0.089379     9 8 0.069477     9 9 -0.224     9 10 0.35521     10 1 -0.89194     10 2 0.13515     10 3 0.6474     10 4 0.88063     10 5 -0.68164     10 6 -0.075601     10 7 -0.63008     10 8 0.168     10 9 -0.84589     10 10 0.46715     1 1 0.98109     1 2 1.5317     1 3 0.46533     1 4 0.12393     1 5 -1.5359     1 6 -0.68207     1 7 -1.6462     1 8 -1.4625     1 9 -1.9343     1 10 -2.3259     2 1 -0.30806     2 2 0.15317     2 3 -0.12825     2 4 0.43423     2 5 -0.71407     2 6 -0.11921     2 7 0.019226     2 8 0.35049     2 9 -1.4383     2 10 -1.5044     3 1 -0.12857     3 2 0.488     3 3 0.53606     3 4 1.0531     3 5 -0.29948     3 6 0.23948     3 7 -0.13268     3 8 -0.24308     3 9 -0.95118     3 10 -0.78948     4 1 -0.11342     4 2 -1.0059     4 3 -0.42124     4 4 -0.34023     4 5 0.078659     4 6 0.88902     4 7 0.23341     4 8 -0.15338     4 9 -0.59989     4 10 -0.89918     5 1 0.16374     5 2 -0.26299     5 3 0.77992     5 4 1.0819     5 5 0.7238     5 6 -0.42438     5 7 -0.19448     5 8 0.08083     5 9 0.48362     5 10 -0.11019     6 1 0.4394     6 2 0.040646     6 3 0.10204     6 4 0.23335     6 5 0.85101     6 6 -0.15969     6 7 -0.11772     6 8 1.4047     6 9 1.591     6 10 0.91408     7 1 0.15331     7 2 -0.61206     7 3 -0.47657     7 4 0.71106     7 5 0.0081887     7 6 -0.040958     7 7 1.0739     7 8 1.175     7 9 1.2231     7 10 0.42095     8 1 0.63768     8 2 0.22313     8 3 0.10201     8 4 0.91221     8 5 -0.2826     8 6 0.53168     8 7 0.043632     8 8 0.12458     8 9 1.1582     8 10 1.036     9 1 -0.83108     9 2 0.31008     9 3 -0.47923     9 4 -0.39486     9 5 0.51966     9 6 0.024501     9 7 0.5772     9 8 -0.18336     9 9 0.51455     9 10 0.99499     10 1 -0.40089     10 2 0.32306     10 3 -0.3433     10 4 0.0047221     10 5 -0.21988     10 6 -0.87257     10 7 -0.22983     10 8 -1.0848     10 9 0.23949     10 10 -0.53886     1 1 -0.89999     1 2 -0.60765     1 3 -1.3965     1 4 -1.1678     1 5 0.50079     1 6 0.7848     1 7 0.31866     1 8 0.33285     1 9 -0.057325     1 10 0.91442     2 1 -1.6385     2 2 -0.66963     2 3 -1.1245     2 4 -0.14442     2 5 0.98717     2 6 0.17885     2 7 0.077735     2 8 -0.22889     2 9 -0.38144     2 10 -0.07872     3 1 -0.77031     3 2 -0.24496     3 3 -0.76457     3 4 -0.74676     3 5 -0.50832     3 6 0.78217     3 7 0.62673     3 8 0.39     3 9 0.065218     3 10 0.488     4 1 0.73499     4 2 -0.22049     4 3 -0.59947     4 4 1.237     4 5 0.95703     4 6 -0.024532     4 7 0.86404     4 8 0.44883     4 9 0.1068     4 10 0.35538     5 1 -0.39373     5 2 0.72013     5 3 0.95173     5 4 0.5653     5 5 1.0025     5 6 -0.1623     5 7 -0.33946     5 8 0.77085     5 9 -0.34686     5 10 0.84909     6 1 0.45291     6 2 1.6368     6 3 0.49003     6 4 1.044     6 5 1.1733     6 6 0.96909     6 7 1.17     6 8 0.789     6 9 0.95848     6 10 0.6586     7 1 0.92626     7 2 1.336     7 3 1.0835     7 4 0.27797     7 5 0.96565     7 6 0.67386     7 7 0.13545     7 8 0.33533     7 9 1.9113     7 10 0.72334     8 1 0.62559     8 2 0.95157     8 3 1.2696     8 4 1.01     8 5 0.75892     8 6 1.1327     8 7 1.0734     8 8 1.5875     8 9 0.27577     8 10 0.52487     9 1 1.3763     9 2 1.1715     9 3 0.33185     9 4 0.7078     9 5 1.7302     9 6 1.6211     9 7 1.2113     9 8 0.28632     9 9 1.4362     9 10 0.05508     10 1 1.6734     10 2 2.0098     10 3 1.3373     10 4 1.9932     10 5 1.6035     10 6 2.3012     10 7 1.4291     10 8 1.1221     10 9 0.66535     10 10 0.92711     1 1 -0.99997     1 2 0.022719     1 3 -0.31203     1 4 0.18867     1 5 0.92715     1 6 1.1044     1 7 0.76145     1 8 -0.17697     1 9 0.25073     1 10 0.31083     2 1 -0.66634     2 2 -1.0159     2 3 -0.72899     2 4 0.32125     2 5 0.54591     2 6 1.1202     2 7 0.64981     2 8 0.03742     2 9 0.77014     2 10 0.012899     3 1 -1.2878     3 2 -1.0884     3 3 -1.0713     3 4 0.22714     3 5 0.36626     3 6 0.20953     3 7 0.14857     3 8 1.3116     3 9 1.2593     3 10 0.65079     4 1 -0.80787     4 2 -0.048305     4 3 -0.87667     4 4 -0.24298     4 5 0.055799     4 6 0.46683     4 7 -0.48911     4 8 -1.1424     4 9 -0.38136     4 10 -1.1061     5 1 -0.27704     5 2 -0.13177     5 3 -0.53845     5 4 -1.1699     5 5 -0.32735     5 6 -0.42001     5 7 0.19617     5 8 0.42689     5 9 -0.13028     5 10 -0.19488     6 1 -0.38319     6 2 0.1336     6 3 -0.49355     6 4 0.14656     6 5 -0.80513     6 6 0.65367     6 7 0.2069     6 8 0.13485     6 9 -0.26875     6 10 0.34008     7 1 -0.62224     7 2 -0.71156     7 3 0.40706     7 4 -0.51364     7 5 -0.045004     7 6 -1.1103     7 7 -1.0858     7 8 0.11086     7 9 -0.97742     7 10 -0.014179     8 1 -1.631     8 2 0.25644     8 3 -0.07203     8 4 0.019954     8 5 -0.78594     8 6 -1.2284     8 7 -0.59174     8 8 -0.29216     8 9 -0.087709     8 10 0.033502     9 1 -1.2977     9 2 -0.18927     9 3 0.67041     9 4 -1.6209     9 5 -0.28613     9 6 -0.025997     9 7 -0.30618     9 8 -1.2943     9 9 -1.587     9 10 -1.9803     10 1 0.46204     10 2 0.074855     10 3 -0.2965     10 4 -0.14154     10 5 -1.1861     10 6 -0.47409     10 7 0.46038     10 8 -0.34689     10 9 -0.80103     10 10 -1.8803     1 1 1.3411     1 2 1.4873     1 3 0.94437     1 4 1.7851     1 5 0.080725     1 6 0.41006     1 7 0.41028     1 8 -0.12254     1 9 0.26331     1 10 -0.42256     2 1 1.2357     2 2 0.32614     2 3 -0.163     2 4 0.98255     2 5 0.72907     2 6 0.31425     2 7 0.45023     2 8 0.26328     2 9 0.029089     2 10 0.95292     3 1 1.073     3 2 1.3613     3 3 0.1275     3 4 0.2476     3 5 0.2768     3 6 0.10945     3 7 0.097306     3 8 -0.39726     3 9 0.33722     3 10 0.061228     4 1 0.48769     4 2 0.78546     4 3 0.84834     4 4 0.97745     4 5 0.42715     4 6 0.077498     4 7 0.01148     4 8 -0.37943     4 9 -1.2956     4 10 0.23022     5 1 1.1574     5 2 0.71203     5 3 1.0309     5 4 1.9035     5 5 0.66869     5 6 0.46976     5 7 -0.31242     5 8 -1.0098     5 9 -0.16017     5 10 -0.055089     6 1 -1.4935     6 2 0.0068701     6 3 0.76397     6 4 0.13412     6 5 -0.36506     6 6 -0.093885     6 7 -0.46345     6 8 0.45901     6 9 -0.3404     6 10 -0.45926     7 1 -0.86785     7 2 -0.1634     7 3 0.043536     7 4 0.62364     7 5 -0.38783     7 6 1.2379     7 7 0.36635     7 8 -0.27407     7 9 0.53955     7 10 0.26264     8 1 -0.49472     8 2 0.33707     8 3 0.34908     8 4 0.39176     8 5 0.49744     8 6 -0.024764     8 7 0.13152     8 8 0.83401     8 9 0.3669     8 10 0.65591     9 1 0.26119     9 2 0.96254     9 3 0.53369     9 4 -0.38361     9 5 -0.31952     9 6 -0.032861     9 7 0.42385     9 8 0.60828     9 9 0.94206     9 10 0.65124     10 1 0.58401     10 2 0.28137     10 3 0.75139     10 4 0.66865     10 5 -0.61473     10 6 0.14     10 7 -0.44916     10 8 -0.41535     10 9 0.072928     10 10 -0.34029     1 1 -1.5582     1 2 -1.9252     1 3 -1.3097     1 4 -2.6259     1 5 -0.86581     1 6 -1.8963     1 7 0.58395     1 8 0.48188     1 9 2.882     1 10 3.7548     2 1 -1.2683     2 2 -1.3114     2 3 -1.1344     2 4 -1.2372     2 5 -1.4385     2 6 -0.38154     2 7 0.21844     2 8 1.2685     2 9 2.0933     2 10 1.6789     3 1 -0.67296     3 2 -0.46436     3 3 -1.1013     3 4 -0.75914     3 5 -0.49573     3 6 -0.071953     3 7 -0.36423     3 8 0.2675     3 9 1.2086     3 10 1.4457     4 1 0.59226     4 2 0.35089     4 3 1.0304     4 4 0.91696     4 5 0.93319     4 6 0.2622     4 7 -0.46816     4 8 0.87504     4 9 0.73656     4 10 0.89449     5 1 0.74513     5 2 0.73003     5 3 1.0218     5 4 0.74639     5 5 0.70416     5 6 -0.29666     5 7 -0.74854     5 8 0.78541     5 9 0.42556     5 10 1.0588     6 1 0.019868     6 2 0.17096     6 3 0.31712     6 4 0.37095     6 5 0.51432     6 6 -0.67542     6 7 -1.0559     6 8 -0.27615     6 9 0.63227     6 10 1.8216     7 1 0.11227     7 2 1.2998     7 3 1.361     7 4 0.84188     7 5 -0.20492     7 6 -0.32484     7 7 -0.58076     7 8 -0.16101     7 9 0.80934     7 10 0.35618     8 1 -0.59692     8 2 1.0501     8 3 0.42325     8 4 -0.78226     8 5 -0.0093265     8 6 0.49451     8 7 -0.57139     8 8 -0.71568     8 9 -1.0949     8 10 -1.322     9 1 -0.3825     9 2 -0.43926     9 3 -0.45886     9 4 0.18823     9 5 0.36943     9 6 -0.01429     9 7 -0.9238     9 8 -0.81771     9 9 -0.58659     9 10 -0.56081     10 1 -1.4616     10 2 -0.83716     10 3 0.83614     10 4 0.28808     10 5 -0.12755     10 6 -1.0188     10 7 -1.2117     10 8 -0.21896     10 9 -1.3982     10 10 -1.9128     1 1 1.5017     1 2 1.1149     1 3 1.5846     1 4 0.49772     1 5 1.6667     1 6 -0.1268     1 7 0.12223     1 8 -0.17968     1 9 0.41108     1 10 1.5414     2 1 0.46015     2 2 -0.18267     2 3 0.18842     2 4 -0.073743     2 5 0.89561     2 6 -0.18944     2 7 0.094503     2 8 0.11414     2 9 0.75025     2 10 0.64978     3 1 0.020306     3 2 -0.25737     3 3 -1.135     3 4 -0.35354     3 5 -0.001421     3 6 -0.12285     3 7 -0.61034     3 8 -0.25232     3 9 0.005621     3 10 -1.3972     4 1 -0.76662     4 2 -0.49325     4 3 -0.76528     4 4 -2.3272     4 5 -0.32207     4 6 0.08859     4 7 0.31069     4 8 0.75689     4 9 -0.36381     4 10 0.10409     5 1 -0.80493     5 2 -0.90191     5 3 -0.64442     5 4 -1.4952     5 5 -0.32429     5 6 -0.017915     5 7 -0.30078     5 8 -0.090466     5 9 1.2656     5 10 0.74413     6 1 0.16902     6 2 0.052467     6 3 -0.23645     6 4 -0.19061     6 5 0.15529     6 6 -0.57483     6 7 -1.1003     6 8 -0.47875     6 9 0.17602     6 10 -0.70837     7 1 0.092632     7 2 -0.73462     7 3 -0.3454     7 4 -0.82274     7 5 -0.51007     7 6 -0.51348     7 7 -1.7752     7 8 -0.58945     7 9 0.30645     7 10 -0.50629     8 1 -0.15482     8 2 -0.4444     8 3 0.1943     8 4 -0.94296     8 5 -0.67503     8 6 -0.7453     8 7 -1.4938     8 8 -0.80872     8 9 -0.052597     8 10 -1.0797     9 1 -0.76072     9 2 0.034279     9 3 0.035288     9 4 -0.47851     9 5 -0.9365     9 6 -0.56425     9 7 -0.75692     9 8 -0.72401     9 9 -0.52698     9 10 -1.2094     10 1 -1.2616     10 2 -0.85022     10 3 -1.2871     10 4 0.61099     10 5 0.058311     10 6 -0.67977     10 7 -0.56729     10 8 -0.97489     10 9 -0.74104     10 10 -0.85617     1 1 -0.81706     1 2 -0.65906     1 3 -0.20288     1 4 0.0097926     1 5 0.74292     1 6 0.31406     1 7 -0.16505     1 8 -0.93741     1 9 -0.22614     1 10 -0.92806     2 1 -0.55144     2 2 -0.059045     2 3 0.095708     2 4 -0.32775     2 5 -0.05043     2 6 -0.90495     2 7 -0.56643     2 8 -0.44872     2 9 -0.5734     2 10 -1.2869     3 1 -0.32205     3 2 -0.49441     3 3 -0.67451     3 4 0.22427     3 5 -0.63732     3 6 -1.7795     3 7 -1.3172     3 8 0.058403     3 9 -0.87751     3 10 -1.176     4 1 0.11018     4 2 -0.67382     4 3 -0.39696     4 4 0.1598     4 5 0.11349     4 6 -0.17633     4 7 -1.0084     4 8 -0.69774     4 9 -0.84491     4 10 -1.569     5 1 0.28874     5 2 -0.28612     5 3 -0.91938     5 4 0.0042778     5 5 0.65028     5 6 0.20493     5 7 -0.23058     5 8 -0.92533     5 9 -1.3078     5 10 -0.050186     6 1 0.82354     6 2 -0.014977     6 3 -0.73324     6 4 0.39265     6 5 -0.10154     6 6 0.43223     6 7 -0.44831     6 8 -0.31979     6 9 -0.37283     6 10 0.035346     7 1 1.1923     7 2 0.95413     7 3 0.98179     7 4 0.39208     7 5 -0.28378     7 6 -1.3313     7 7 -1.5069     7 8 -1.0779     7 9 -0.62778     7 10 1.8462     8 1 0.68447     8 2 1.1256     8 3 -0.1701     8 4 0.85436     8 5 -0.0047328     8 6 -0.66924     8 7 -0.34008     8 8 -0.43798     8 9 -0.47153     8 10 -0.11979     9 1 1.1134     9 2 1.5582     9 3 0.49476     9 4 0.46428     9 5 0.9354     9 6 0.84867     9 7 0.75609     9 8 -1.0049     9 9 -0.45265     9 10 -0.74576     10 1 -0.27515     10 2 0.35115     10 3 1.1861     10 4 0.98898     10 5 0.21486     10 6 -0.83618     10 7 0.015597     10 8 -0.83299     10 9 -0.2658     10 10 0.061564     1 1 1.0847     1 2 1.9526     1 3 1.4559     1 4 1.2725     1 5 0.11601     1 6 -0.55716     1 7 -1.1344     1 8 -0.92588     1 9 0.14647     1 10 -0.54407     2 1 1.8636     2 2 1.044     2 3 -0.053585     2 4 0.32422     2 5 -0.17135     2 6 -0.16723     2 7 -1.0589     2 8 -0.2573     2 9 -0.2436     2 10 -0.86413     3 1 1.4981     3 2 0.66796     3 3 -0.010454     3 4 0.016212     3 5 -0.23302     3 6 -0.59975     3 7 -1.7644     3 8 -0.51085     3 9 -0.3481     3 10 -2.3364     4 1 1.2338     4 2 0.59601     4 3 0.10002     4 4 -0.57298     4 5 0.27338     4 6 -0.63656     4 7 -0.70863     4 8 -1.3946     4 9 -0.85695     4 10 -1.4846     5 1 1.0698     5 2 0.397     5 3 1.1341     5 4 0.90212     5 5 0.52142     5 6 -0.24757     5 7 -0.71933     5 8 -1.2248     5 9 -1.0037     5 10 -1.4199     6 1 0.697     6 2 0.67897     6 3 -0.037977     6 4 0.42927     6 5 0.73508     6 6 0.46973     6 7 -0.77912     6 8 -1.5954     6 9 0.078578     6 10 -1.1593     7 1 -0.88098     7 2 0.35062     7 3 0.16084     7 4 0.41738     7 5 0.74172     7 6 0.46485     7 7 0.075388     7 8 -0.84538     7 9 -0.53554     7 10 -1.5093     8 1 0.58167     8 2 0.16468     8 3 0.34762     8 4 0.16298     8 5 -0.27006     8 6 0.08014     8 7 1.5419     8 8 0.27913     8 9 -0.69523     8 10 -0.48911     9 1 0.7055     9 2 0.078586     9 3 -0.28458     9 4 -0.37242     9 5 -0.21641     9 6 -0.34865     9 7 0.85598     9 8 0.18734     9 9 0.22394     9 10 1.0592     10 1 -0.092963     10 2 1.3653     10 3 1.2755     10 4 -0.14577     10 5 -0.36608     10 6 -0.0095181     10 7 0.47785     10 8 0.49047     10 9 0.52567     10 10 1.3327     1 1 0.084506     1 2 -0.13379     1 3 -0.47491     1 4 -0.043682     1 5 -0.65521     1 6 -0.052312     1 7 -0.12092     1 8 -1.0912     1 9 -0.66403     1 10 0.08132     2 1 -0.045233     2 2 0.48315     2 3 0.096673     2 4 -0.33502     2 5 -0.52595     2 6 -0.68735     2 7 -0.34867     2 8 -1.2754     2 9 -1.3937     2 10 -0.4887     3 1 -0.31061     3 2 -0.022152     3 3 0.47162     3 4 0.20441     3 5 0.059868     3 6 -0.2901     3 7 -0.4409     3 8 -1.1976     3 9 -0.73539     3 10 -1.435     4 1 0.30343     4 2 0.11239     4 3 1.8881     4 4 1.2312     4 5 0.45089     4 6 0.084228     4 7 -0.19204     4 8 -1.1114     4 9 -1.7262     4 10 -0.56662     5 1 -0.17351     5 2 0.34411     5 3 0.85914     5 4 0.43333     5 5 0.11873     5 6 0.41948     5 7 0.49702     5 8 -0.67355     5 9 -0.71335     5 10 -0.56962     6 1 -0.72632     6 2 0.182     6 3 0.29566     6 4 0.21552     6 5 -0.46146     6 6 -0.059419     6 7 0.11557     6 8 -1.0477     6 9 -0.29064     6 10 -1.2339     7 1 -0.27475     7 2 1.0619     7 3 0.40953     7 4 -0.050585     7 5 0.00087871     7 6 0.40448     7 7 -0.25119     7 8 -0.41342     7 9 -0.133     7 10 -0.28252     8 1 0.23176     8 2 2.0917     8 3 1.5435     8 4 0.075721     8 5 -0.10405     8 6 -0.063779     8 7 1.1413     8 8 0.40685     8 9 -0.58853     8 10 -0.49277     9 1 0.04969     9 2 1.481     9 3 0.69981     9 4 0.14286     9 5 0.065013     9 6 -0.37731     9 7 0.40477     9 8 0.68332     9 9 -0.18014     9 10 -1.6154     10 1 0.30045     10 2 -0.43197     10 3 0.20984     10 4 0.40885     10 5 -0.012176     10 6 0.012924     10 7 -0.97699     10 8 0.25918     10 9 -1.2407     10 10 -2.4877     1 1 1.9568     1 2 1.4     1 3 0.43733     1 4 -0.43196     1 5 -0.16741     1 6 0.1497     1 7 0.21556     1 8 1.1755     1 9 1.531     1 10 1.9811     2 1 1.0004     2 2 0.50312     2 3 0.57932     2 4 0.34859     2 5 0.35882     2 6 0.5033     2 7 -0.30434     2 8 0.53376     2 9 0.88751     2 10 0.82788     3 1 0.28532     3 2 -0.68586     3 3 0.13507     3 4 0.0048657     3 5 -0.12591     3 6 -0.32573     3 7 0.17317     3 8 0.67218     3 9 -0.68763     3 10 0.35134     4 1 -0.31966     4 2 -0.51346     4 3 0.26837     4 4 -0.037354     4 5 0.45119     4 6 -0.57281     4 7 -0.11622     4 8 0.076975     4 9 0.60791     4 10 0.66343     5 1 0.18332     5 2 0.2538     5 3 0.34382     5 4 0.5546     5 5 0.21333     5 6 0.80252     5 7 0.47662     5 8 1.1387     5 9 0.39022     5 10 0.76584     6 1 -0.63939     6 2 0.47384     6 3 -0.67462     6 4 0.3125     6 5 0.70997     6 6 0.59881     6 7 0.99534     6 8 1.1245     6 9 0.6379     6 10 0.3624     7 1 0.11995     7 2 -0.78728     7 3 -0.13309     7 4 -0.38062     7 5 0.43758     7 6 0.029345     7 7 0.40294     7 8 -0.32999     7 9 0.83065     7 10 0.85057     8 1 1.3007     8 2 0.66175     8 3 1.1035     8 4 0.26512     8 5 0.15     8 6 0.67162     8 7 0.60928     8 8 1.1583     8 9 0.86711     8 10 0.052528     9 1 0.72511     9 2 0.83225     9 3 1.4767     9 4 1.2465     9 5 0.62551     9 6 0.18133     9 7 0.54219     9 8 0.61455     9 9 0.069688     9 10 0.91685     10 1 0.87253     10 2 2.1886     10 3 1.4976     10 4 1.3471     10 5 1.614     10 6 0.81966     10 7 0.4422     10 8 -0.011136     10 9 -0.6524     10 10 -0.39135     1 1 1.0954     1 2 0.68865     1 3 0.31133     1 4 -0.47359     1 5 0.16099     1 6 1.3591     1 7 -0.32577     1 8 -0.13171     1 9 0.62947     1 10 0.44251     2 1 0.92282     2 2 0.91816     2 3 0.76012     2 4 0.24178     2 5 -0.068355     2 6 0.25886     2 7 0.22525     2 8 -0.79252     2 9 0.14529     2 10 1.1726     3 1 0.38776     3 2 0.068483     3 3 0.058235     3 4 0.24309     3 5 -0.45677     3 6 0.51896     3 7 1.3932     3 8 0.64736     3 9 0.48068     3 10 0.74987     4 1 0.30978     4 2 -0.098167     4 3 -0.10243     4 4 -0.060663     4 5 0.66854     4 6 0.11733     4 7 1.2778     4 8 0.21566     4 9 0.46415     4 10 0.63132     5 1 0.22002     5 2 0.27927     5 3 0.39181     5 4 1.0721     5 5 0.8424     5 6 -0.75326     5 7 -0.1773     5 8 0.31369     5 9 0.78829     5 10 0.86142     6 1 0.60782     6 2 1.1523     6 3 -0.24421     6 4 -0.70425     6 5 -0.91406     6 6 -0.25592     6 7 -0.31783     6 8 0.90311     6 9 1.0813     6 10 1.6632     7 1 0.23327     7 2 0.13924     7 3 0.24131     7 4 -0.75094     7 5 0.16098     7 6 0.58843     7 7 0.14414     7 8 1.1983     7 9 1.2397     7 10 1.282     8 1 0.054621     8 2 0.96913     8 3 1.1722     8 4 0.88584     8 5 0.43832     8 6 0.35612     8 7 0.14986     8 8 0.2759     8 9 0.561     8 10 0.1761     9 1 0.47944     9 2 1.237     9 3 1.5481     9 4 0.66968     9 5 -0.3433     9 6 0.0070064     9 7 -0.66577     9 8 0.091003     9 9 0.47818     9 10 0.062905     10 1 1.0756     10 2 0.84519     10 3 1.2029     10 4 0.63381     10 5 -0.20783     10 6 0.23517     10 7 0.29097     10 8 -0.29507     10 9 0.0818     10 10 -0.061862     1 1 -0.60204     1 2 0.65     1 3 -0.33662     1 4 -0.26463     1 5 0.65212     1 6 -0.6876     1 7 -0.94691     1 8 -1.0972     1 9 -0.19294     1 10 -0.28093     2 1 -0.87526     2 2 -0.081139     2 3 -0.87035     2 4 -0.74718     2 5 -0.26173     2 6 -1.1377     2 7 -0.81349     2 8 0.013986     2 9 -0.56303     2 10 0.033915     3 1 -0.80039     3 2 -0.38038     3 3 -0.14681     3 4 -1.3574     3 5 -1.1031     3 6 -0.53416     3 7 -0.71843     3 8 -0.30968     3 9 -0.17786     3 10 -0.32835     4 1 0.073412     4 2 -0.21385     4 3 -0.2436     4 4 -0.59427     4 5 0.10118     4 6 -0.92508     4 7 0.32233     4 8 -0.038958     4 9 0.074636     4 10 -0.7109     5 1 0.626     5 2 -0.48074     5 3 -0.71395     5 4 -1.6909     5 5 -0.56445     5 6 -0.83992     5 7 0.0088717     5 8 -0.34586     5 9 -0.17763     5 10 0.85025     6 1 0.72631     6 2 -0.83441     6 3 -0.85932     6 4 -1.2019     6 5 -0.3461     6 6 -0.52261     6 7 -0.2188     6 8 0.1562     6 9 0.1176     6 10 -0.0045701     7 1 0.90036     7 2 -0.66526     7 3 -1.5311     7 4 -0.70255     7 5 -0.39794     7 6 0.20198     7 7 0.45935     7 8 0.4518     7 9 0.82838     7 10 0.76867     8 1 0.16852     8 2 -0.63111     8 3 -1.1804     8 4 -0.69962     8 5 0.05342     8 6 0.16314     8 7 0.59872     8 8 -0.12264     8 9 -0.73362     8 10 0.057699     9 1 0.64873     9 2 -0.052743     9 3 -0.5923     9 4 -0.2954     9 5 -0.74368     9 6 0.32805     9 7 0.079574     9 8 -0.63894     9 9 -0.21909     9 10 0.64924     10 1 0.17714     10 2 -1.1491     10 3 -0.45945     10 4 -0.85382     10 5 1.4735     10 6 -0.029148     10 7 -0.4445     10 8 -0.72929     10 9 -0.10109     10 10 -0.027009


title: Interpreting a simple linear regression coefficient scaled by the mean of y?

text: I have come across an industry example of a simple linear regression ($y=a+bx+\epsilon$) where the slope coefficient has been adjusted by the mean of $y$ ($b/\text{mean}(y)$) and described as a "slope impact". The exact interpretation that they use is that this slope impact represents "the percentage change in $y$ for every unit above the average value of $x$". I am not sure this interpretation is correct but would appreciate any feedback you could give on this calculation.


title: Statistical methods for analysing correlation between stock price index and natural disasters

text: Can someone help me identify what statistical method (or any method) that I can use to correlate the effects of natural disasters on Stock Market Index?


title: Threshold for psychometric experiment with three (or more) possible choices

text: I'm running an experiment where subjects are viewing a stimulus, and they have to decide if the stimulus is brighter/darker/identical in respect to a standard image. Stimulus luminosity actually varies, and it can be brighter, darker or identical to the standard. What i' triyng to estimate is the threshold (or the just noticeable difference) needed to correctly identify (let's say 75% of times) the stimulus as darker or brighter then the standard. I prefer to run the experiment with the method of constant stimuli, randomizing the standard luminosity. From dummy data, I suppose that the plot of probability of response will follow the graph here showed. On Y-axis the probability; on X-axis the luminosity of the stimulus (-7 is darker, +7 is brighter, 0 is identical luminosity). The red curve identify the probability of answering "darker"; the green curve the probability of answering "brighter", and the blue line the probability of answering "identical luminosity". ![enter image description here](http://i.stack.imgur.com/gYzvI.jpg) How can I consider an experiment like this? What is the correct model of analysis? From google-searches I've found that it could be considered a "multinomial ordered logit model", but I'm not sure. And by the way, I was not able to find how to determinate the threshold or the JND in experiments like that. Any references are welcome.


title: Dealing with multicollinearity in logistic regresson model

text: I have 5 quantitative predictor variables in my logit model, and when I use the `cor` function in R on those 5 variables, I see that $x_{1}$ and $x_{2}$ have correlation coefficient $r_{12}=0.718$ and $x_{3}$ and $x_{4}$ have $r_{34}=0.654$. The VIF values are all less than or equal to 2. So does this suggest that there's no multicollinearity in my model? If there is, my question is how can I deal with multicollinearity? Should I create interaction terms? Or would factor analysis/principal component analysis be applicable here? I'm not familiar with those concepts but I can definitely read up on them.


title: How do I find similarities between two sets of data

text: I have a group of data with 12 different football players, and they are rated for 11 different skills (speed, skill, flair, etc). I am looking to pair up individuals based on similar footballers, and was thinking of just taking the average for each player and matching them up. Does anyone have any better suggestions, i.e. a statistical test to see which players are best matched. And from this we can have fair teams... Hope this makes sense,               +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     |             |player|player|player|player|player|player|player|player|player|player|player|player|     |             |     1|     2|     3|     4|     5|     6|     7|     8|     9|    10|    11|    12|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Speed       |  12.1|  15.9|  10.9|  12.4|   7.6|  13.6|   9.3|   6.9|  13.4|   8.6|   4.7|  13.1|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Skill       |   8.1|  14.1|  13.6|   7.7|   6.0|  10.4|   8.1|   7.0|  13.1|   9.7|  17.4|  14.7|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Flair       |   6.6|  13.3|  11.3|   7.4|   5.9|  10.0|   8.1|   8.4|  13.6|   8.1|  16.6|  15.0|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Finishing   |   9.4|  13.3|  11.3|   7.4|   7.4| 12.00|  10.6|  11.9|  13.4|  17.4|  13.6|  12.3|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Long  Shots |  10.1|  13.1|  11.6|   7.4|   6.6|   8.7|  10.7|  13.6|  13.4|  15.6|  12.4|   9.6|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Crossing    |  10.0|  12.7|  15.7|   9.6|   6.6|  10.3|  11.3|   9.9|  11.7|   8.9|  13.6|  11.1|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Dribbling   |  12.0|  15.0|  13.3|   9.3|   7.1|  11.0|  10.6|   8.4|  12.7|  10.3|  16.1|  13.9|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Tackling    |  18.4|  14.7|   8.3|  15.1|   9.3|   5.0|  14.6|  11.9|  13.0|   7.4|   9.6|   5.9|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Heading     |  14.7|  15.9|   8.4|  11.6|   7.7|   5.0|  13.1|   9.6|  13.1|  13.1|   8.0|   6.3|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     |Determination|  19.0|  14.9|   9.1|  15.6|  11.7|   7.9|  14.4|  10.9|  12.9|  12.1|  11.4|   5.1|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Passing     |  12.9|  15.1|  15.9|   9.6|   9.3|  11.3|  12.1|  10.9|  16.3|  11.1|  17.6|  13.3|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Total       | 121.3| 142.1| 118.4| 100.7|  77.6|  91.6| 113.7| 102.3| 133.3| 113.9| 136.3| 107.1|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Maximum     |   220|   220|   220|   220|   220|   220|   220|   220|   220|   220|   220|   220|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+     | Index rating| 55.13| 64.61| 53.83| 45.78| 35.26| 41.62| 51.69| 46.49| 60.58| 51.75| 61.95| 48.70|     +-------------+------+------+------+------+------+------+------+------+------+------+------+------+


title: Obtaining sequence of lambda values for training glmnet model via `caret`

text: I have multiple models that I'm training using `train` in the `caret` package, all while using the same cross validation folds to compare their performance. I'd like to use glmnet via caret's `train` and specify the sequence of lambda values by obtaining them like how the `cv.glmnet`function inside the `glmnet` package does it. The documentation for `glmnet` states that for `cv.glmnet`: > The function runs glmnet nfolds+1 times; **the first to get the lambda > sequence** , and then the remainder to compute the fit with each of the > folds omitted... Here's a prototype of what the situation using `caret` looks like:               ##Create cross-validation folds to use for multiple models     set.seed(9)     MyFolds <- createMultiFolds(trainset$Response, times = 100)           ##my model train control function     fitControl1 <- trainControl(method = "repeatedcv",                                 repeats = 100,                                 index = MyFolds,                                 classProbs = TRUE,                                 allowParallel = TRUE,                                 summaryFunction = fiveStats)          ## fit the glmnet model     model <- train(FORMULA, data=trainset,                    method='glmnet',                    metric = "ROC",                    tuneGrid = expand.grid(.alpha=c(0,1),.lambda=seq(0,.25,by=0.005)),                    trControl=fitControl1)     model      Inside `tuneGrid = expand.grid(.alpha=c(0,1),.lambda=seq(0,.25,by=0.005))`, I'd like to obtain a sequence of lambda values similar to how `cv.glmnet` does it, but using the cv-folds I've already specified via:               set.seed(9)     MyFolds <- createMultiFolds(trainset$Response, times = 100)      How do I achieve this?


title: Alternatives to one-way ANOVA for heteroskedastic data

text: I have data from 3 groups of algae biomass ($A$, $B$, $C$) which contain unequal sample sizes ($n_A=15$, $n_B=13$, $n_C=12$) and I would like compare if these groups are from the same population. One-way ANOVA would definitely be the way to go, however upon conducting normality tests on my data, heteroskedascity seems to the main issue. My raw data, without any transformation, produced a ratio of variances ($F_{\max} = 19.1$) which is very much higher than the critical value ($F_{\rm crit} = 4.16$) and therefore I cannot perform one-way ANOVA. I also tried transformation to normalise my data. Even after trials of various transformations (log, square root, square), the lowest $F_{\max}$ produced after transformation with a $\log_{10}$ transformation was $7.16$, which was still higher compared to $F_{\rm crit}$. Can anyone here advise me on where to go from here? I can't think of other methods of transformation to normalise by data. Are there any alternatives to a one-way ANOVA? P.S.: my raw data are below:               A: 0.178 0.195 0.225 0.294 0.315 0.341 0.36  0.363 0.371 0.398 0.407 0.409 0.432         0.494 0.719     B: 0.11  0.111 0.204 0.416 0.417 0.441 0.492 0.965 1.113 1.19  1.233 1.505 1.897     C: 0.106 0.114 0.143 0.435 0.448 0.51  0.576 0.588 0.608 0.64  0.658 0.788 0.958


title: Estimating parameters of an unknown PID controller

text: Say that I have your standard PID controller at work. To keep it extremely simple imagine I have a target $x^*$ on the variable $x$. Then the controller is: $y(t) = K_p ( x^* - x_t) + K_i \int_0^t (x^* - x_{\tau})d \tau + K_d \frac{d}{dt} (x^* -x_t)$ Now, I do have numerical observations for all $t \in [0,T]$ of both $y_t$ and $x_t$. This data generating process is a PID controller of which I do not know the three parameters. I'd like to estimate them! Now if the data was perfect I would just need the first 3 observations to solve a system of equation and I'd be done. The data isn't perfect though (imagine an $\epsilon$ added to my observed $y_t$). How should I proceed? I couldn't find literature regarding this problem so I am somewhat at a loss. How do you attack this?   * Do I perform an OLS? I can approximate the integrative and derivative part with my data set but I wonder how the errors will compound. Plus there is $t$ and I wonder if any time series object would do the job better   * I could try brute search or some weird metropolis hastings trick where I keep trying random parameters and choose the ones that minimize SSR or prediction error, but I think that's my last option. Any ideas?


title: STL on time series with missing values for anomaly detection

text: I am trying to detect anomalous values in a time series of climatic data with some missing observations. Searching the web I found many available approaches. Of those, stl decomposition seems appealing, in the sense of removing trend and seasonal components and studying the remainder. Reading http://cs.wellesley.edu/~cs315/Papers/stl%20statistical%20model.pdf , `stl` appears to be flexible in determining the settings for assigning variability, unaffected by outliers and possible to apply despite missing values. However, trying to apply it in `R`, with four years of observations and defining all the parameters according to http://stat.ethz.ch/R-manual/R-patched/library/stats/html/stl.html , I encounter error: `"time series contains internal NAs"` (when `na.action=na.omit`), and   `"series is not periodic or has less than two periods"` (when `na.action=na.exclude`). I have double checked that the frequency is correctly defined. I have seen relevant questions in blogs, but didn't find any suggestion that could solve this. Is it not possible to apply `stl` in a series with missing values? I am very reluctant to interpolate them, as I do not want to be introducing (and consequently detecting...) artifacts. For the same reason, I do not know how advisable it would be to use ARIMA approaches instead (and if missing values would still be a problem). Please share if you know a way to apply `stl` in a series with missing values, or if you believe my choices are methodologically not sound, or if you have any better suggestion. I am quite new in the field and overwhelmed by the heaps of (seemingly...) relevant information.


title: Minimum population for decision trees and association rules

text: Hi I'm quite new to this and I'm playing around with R and Microsoft's SSAS. Does anyone have a rule of thumb how big a data set has to be for association rules and decision trees to be statistically valid. In uni our assignments use data sets with about a 1000 in.


title: Clustered standard errors in 2-period Dif-in-Dif?

text: in order to rectify invalid t-stats because of autocorrelation in Difference- in-Differences (DnD) models, Duflo et al (2004) propose (among other solutions) to collapse data so as to have a before-after DnD. My question: What if we already HAVE a before-after setting (e.g. pre-reform vs post- reform)? Would it be still ok to use clustered SE (provided the number of clusters would be high enough)? Or is it automatically wrong to use them EVEN IF we already have data that has 2 time periods only. Would appreciate any comments. Thanks.


title: I need to check if I'm correct, with conditional and bayes

text: I have this problem:               A Ph.D. graduate has applied for a job with two      universities: A and B. The graduate feels that she has a      60% chance of receiving an offer from university A and a      50% chance of receiving an offer from university B. If      she receives an offer from university B, she believes that      she has an 80% chance of receiving an offer from      university A.     a) What is the probability that both universities will make      her an offer?     b) What is the probability that at least one university will      make her an offer?     c) If she receives an offer from university B, what is the      probability that she will not receive an offer from      university A?      And I have done it like this:                   A         A'     B  0.4        0.1    p(B)=0.5              B' 0.2        0.3   p(B')=0.5       p(A)=0.6   p(A')=0.4          p(A|B)=0.8     p(A and B)=p(B)p(A|B)=0.5(0.4)=0.4          a) p(both)=p(A and B)=0.4     b) p(at least 1)=1-p(none)= 1-p(A' and B')=0.7     c)p(A'|B)=p(B|A')p(A')/p(B)=(0.25(0.4))/0.5=0.2          p(B|A')=p(A' and B)/p(A')= 0.1/0.4=0.25      So is this the correct way of approaching it? Is it correct? Thanks.


title: stepwise selection on Negative Binomial regression model

text: I know that I can perform a stepwise selection on ordinary linear regression model based on the t-value. but what about Negative Binomial regression model - or GLM in general? Does it theoretically make sense to perform a stepwise selection on variables bases on P-value resulting from "z-scores" when I am trying to improve Negative Binomial regression model? Is such method recommended? thank you


title: confidence band multiple regression

text: I have realized a multiple linear regression with R. When I realize the graphs diagnostics, I would like to draw the confidence band along with the loess curve in the first plot (residual vs. fitted values). How can I do that with R (R codes) ? Best,


title: Randomized block design for repeated trials

text: Suppose I measure reaction time of 4 students, initially before any treatment, then after treatment 1, and finally after treatment 2. (Note all treatments will eventually be applied to all students.) Also because of potential error with reaction time measurements, suppose 3 trials are done for each student- treatment combination. Now if there was only 1 student-treatment combination, a simple repeated measures ANOVA could be applied. In this case, is it statistically OK to instead use a randomized block design, where each block is comprised of the trials relating to a particular student?


title: Explanation for unstable lasso regression coefficients?

text: I have run a lasso regression on a dataset of 100 observations and 80 variables (using 10-fold cross-validation to find the minimum lambda subsequently used in the final model). The lasso regression found approximately 40 of the variables to have non-zero coefficients. I wanted to check my model and therefore divided those 100 observations into two sets (70/30 - the idea being that I would have a train and test set) and ran a lasso regression on the 70. All coefficients calculated were 0 (except for the intercept) - a dramatically different result than in the first model using all 100 variables. Confused, I ran another lasso regression on just the 30 observations and found 3 variables to have non-zero coefficients. I assume that my drastically different results stem from the fact that the data I have does not do a good job of explaining the dependent variable, but perhaps there is a better explanation? In case this is helpful - I am interested in using lasso regression for prediction.


title: Estimating model error in $k$-nearest neighbours with strongly spatially autocorrelated training data

text: In the palaeoclimate world, palaeoecologists have used spatial training sets of say sea-surface temperture (SST) and related this to micro-organisms living at the locations where SST was measured. A popular model to predict SST from species assemblages (so you can count the same species from sediment core and thus predict SST back in time) is the Modern Analogue Technique (MAT). MAT is our name for $k$ nearest neighbours, but we are using it in a regression rather than the more commonly encountered classification setting. In MAT ($k$-NN) we select the $k$ closest samples in terms of species composition to infer SST. Under standard (leave-one-out, $k$-fold, or bootstrapping) crossvalidation these models have very good performance stats and low model error compared to other modelling approaches. It has been shown that much of this improvement is due to the large degree of spatial autocorrelation in the training data. MAT is more able to use this information than some other techniques, hence the favourable RMSEP). One solution suggested to this problem (of getting a realistic RMSEP or at least a better understanding of model error/uncertainty or even quality) is to use $h$-block crossvalidation, where samples with distance $h$ of the target sample are removed during crossvalidation. However, with MAT, these will necessarily be the observations that the model needs in order to predict that SST; by throwing out the samples within distance $h$ you are excluding sites with similar species composition and hence forcing the model to extrapolate. This is much less of a problem where the relationships between species and SST are "modelled" for example. I'm wondering if there are other ways in which spatial autocorrelation can be handled within a $k$-NN framework. Are there better/alternative ways to tackle the problem of estimating an RMSEP for a $k$-NN model in a situation as I describe (strong spatial dependence in the training data). I have been thinking about using a dumb null model using spatial distance against which to compare the $k$-NN model, but that doesn't get at estimating an RMSEP that is not biased due to the spatial autocorrelation.