2015-09-04 4 views
4

В «Элементы статистического обучения» по Tibshirani, при сравнении методом наименьших квадратов/линейных моделей и Knn эти 2 сценария заявлены:Как генерировать данные для гауссовых распределений в этих двух сценариях в R?

Сценарий 1: Тренировочные данные в каждом классе были получены из двумерных распределений Гаусса с коррелированы компонентов и различных средств.

Сценарий 2. Данные обучения в каждом классе были получены из смеси 10 низкочастотных гауссовских распределений с самими отдельными средствами , распределенными как гауссовские.

Идея заключается в том, что первая лучше подходит для наименьших квадратов/линейных моделей, а второй для Knn как модели (те, с более высокой дисперсией от того, что я понимаю, так как КНН учитывает ближайшие точки и не все точки) ,

В R, как бы я смоделировал данные для обоих сценариев?

Конечная цель состоит в том, чтобы воспроизвести оба сценария, чтобы доказать, что эффективно 1-й лучше объясняет линейная модель, чем вторая.

Спасибо!

+0

Хотя этот вопрос имеет статистический компонент, было бы интересно объяснить различие между два сценария - его акцент на «R» заставил его собирать чисто программно-ориентированные ответы, делая его вне темы на CV. – whuber

+0

Я действительно хочу, чтобы иметь возможность воспроизвести 2 сценария, чтобы лучше понять различия и как они применяются к линейным моделям. – Sofia

ответ

1

Это может быть сценарий 1

library(mvtnorm) 

N1 = 50 
N2 = 50 
K = 2 

mu1 = c(-1,3) 
mu2 = c(2,0) 

cov1 = 0 
v11 = 2 
v12 = 2 
Sigma1 = matrix(c(v11,cov1,cov1,v12),nrow=2) 

cov2 = 0 
v21 = 2 
v22 = 2 
Sigma2 = matrix(c(v21,cov2,cov2,v22),nrow=2) 

x1 = rmvnorm(N1,mu1,Sigma1) 
x2 = rmvnorm(N2,mu2,Sigma2) 

Это может быть кандидатом для моделирования из гауссовой смеси:

BartSimpson <- function(x,n = 100){ 
    means <- as.matrix(sort(rnorm(10))) 
    dens <- .1*rowSums(apply(means,1,dnorm,x=x,sd=.1)) 
    rBartSimpson <- c(apply(means,1,rnorm,n=n/10,sd=.1)) 
    return(list("thedensity" = dens,"draws" = rBartSimpson)) 
} 

x <- seq(-5,5,by=.01) 

plot(x,BartSimpson(x)$thedensity,type="l",lwd=4,col="yellow2",xlim=c(-4,4),ylim=c(0,0.6)) 
+0

Спасибо :) Почему этот вид? – Sofia

+0

вполне может быть лишним, я что-то исправлял и пытался это сделать, но это было не главное, как выяснилось –

1

В приведенном ниже коде я сначала создаю 10 различных средств классов, а затем использую средства для получения случайных значений из этих средств. Код идентичен для двух сценариев, но вам нужно будет отрегулировать дисперсию внутри и между классами, чтобы получить желаемые результаты.

Сценарий 1:

Здесь вы хотите создать 10 классов с различными средствами (я предполагаю, что средства следуют двумерное гауссово распределение). Разница между классами намного меньше, чем разница между классами.

library(MASS) 
n <- 20 
# subjects per class 
classes <- 10 
# number of classes 
mean <- 100 
# mean value for all classes 
var.between <- 25 
# variation between classes 
var.within <- 225 
# variation within classes 
covmatrix1 <- matrix(c(var.between,0,0,var.between), nrow=2) 
# covariance matrix for the classes 
means <- mvrnorm(classes, c(100,100), Sigma=covmatrix1) 
# creates the means for the two variables for each class using variance between classes 
covmatrix2 <- matrix(c(var.within,0,0,var.within), nrow=2) 
# creates a covariance matrix for the subjects 
class <- NULL 
values <- NULL 
for (i in 1:10) { 
    temp <- mvrnorm(n, c(means[i], means[i+classes]), Sigma=covmatrix2) 
    class <- c(class, rep(i, n)) 
values <- c(values, temp) 
} 
# this loop uses generates data for each class based on the class means and variance within classes 
valuematrix <- matrix(values, nrow=(n*classes)) 
data <- data.frame (class, valuematrix) 
plot(data$X1, data$X2) 

В качестве альтернативы, если вы не заботитесь о задании дисперсии между классами, и вы не хотите какой-либо корреляции внутри классов, вы можете просто сделать это:

covmatrix <- matrix(c(225, 0, 0, 225), nrow=2) 
# specifies that the variance in both groups is 225 and no covariance 
values <- matrix(mvrnorm(200, c(100,100), Sigma=covmatrix), nrow=200) 
# creates a matrix of 200 individuals with two values each. 

Сценарий 2:

Здесь единственное отличие состоит в том, что разница между классами больше, чем вариация внутри классов. Попробуйте обмениваясь значение переменной var.between примерно до 500, а переменная var.within до 25, и вы увидите четкую кластеризацию рассеивания:

n <- 20 
# subjects per class 
classes <- 10 
# number of classes 
mean <- 100 
# mean value for all classes 
var.between <- 500 
# variation between classes 
var.within <- 25 
# variation within classes 
covmatrix1 <- matrix(c(var.between,0,0,var.between), nrow=2) 
# covariance matrix for the classes 
means <- mvrnorm(classes, c(100,100), Sigma=covmatrix1) 
# creates the means for the two variables for each class using variance between classes 
covmatrix2 <- matrix(c(var.within,0,0,var.within), nrow=2) 
# creates a covariance matrix for the subjects 
class <- NULL 
values <- NULL 
for (i in 1:10) { 
    temp <- mvrnorm(n, c(means[i], means[i+classes]), Sigma=covmatrix2) 
    class <- c(class, rep(i, n)) 
values <- c(values, temp) 
} 
# this loop uses generates data for each class based on the class means and variance within classes 
valuematrix <- matrix(values, nrow=(n*classes)) 
data <- data.frame (class, valuematrix) 
plot(data$X1, data$X2) 

Участок должен подтвердить, что данные сгруппированы.

Надеюсь, это поможет!

+0

Я ожидал, что сценарий 1 будет лучше объяснен применением линейной модели (model = lm (X2 ~ X1, data = data)), чем сценарий 2, но это не кажется истинным, хотя я вижу, что в сценарии 2 данные кластеризуется. Даже изменение var.betweeen и var.within в сценарии 1 до низких значений (например, 2), похоже, не обеспечивает последовательного предоставления более качественных p-значений для линейной модели в сценарии 1. Знаете ли вы, почему? – Sofia

+0

Извините, но я мало знаю о моделях knn. В основном я пытался помочь в создании данных. – JonB

0

С помощью обоих ответов здесь я в конечном итоге с помощью этого:

mixed_dists = function(n, n_means, var=0.2) { 
    means = rnorm(n_means, mean=1, sd=2) 
    values <- NULL 
    class <- NULL 
    for (i in 1:n_means) { 
     temp <- rnorm(n/n_means, mean=means[i], sd=0.2) 
     class <- c(class, rep(i, n/n_means)) 
     values <- c(values, temp) 
    } 
    return(list(values, class)); 
} 

N = 100 

#Scenario 1: The training data in each class were generated from bivariate Gaussian distributions 
#with uncorrelated components and different means. 
scenario1 = function() { 
    var = 0.5 
    n_groups = 2 
    m = mixed_dists(N, n_groups, var=var) 
    x = m[[1]] 
    group = m[[2]] 
    y = mixed_dists(N, n_groups, var=var)[[1]] 
    data = matrix(c(x,y, group), nrow=N, ncol=3) 
    colnames(data) = c("x", "y", "group") 
    data = data.frame(data) 
    plot(x=data$x,y=data$y, col=data$group) 
    model = lm(y~x, data=data) 
    summary(model) 
} 



#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. 
scenario2 = function() { 
    var = 0.2 # low variance 
    n_groups = 10 
    m = mixed_dists(N, n_groups, var=var) 
    x = m[[1]] 
    group = m[[2]] 
    y = mixed_dists(N, n_groups, var=var)[[1]] 
    data = matrix(c(x,y, group), nrow=N, ncol=3) 
    colnames(data) = c("x", "y", "group") 
    data = data.frame(data) 
    plot(x=data$x,y=data$y, col=data$group) 
    model = lm(y~x, data=data) 
    summary(model) 
} 
# scenario1() 
# scenario2() 

Таким образом, данные в сценарии 1 четко разделены в 2 классах, а данные i n сценарий 2 имеет около 10 кластеров и не может быть четко разделен с использованием прямой линии. На самом деле, используя линейную модель для обоих сценариев, можно видеть, что в среднем она будет лучше применяться к сценарию 1, чем к сценарию 2.

 Смежные вопросы

  • Нет связанных вопросов^_^