Se você desejar entender o significado de correlação experimente:

n=10000000
cor_x1x2=-0.2
medx1=0.0067
varx1=0.0017
sdx1=varx1^0.5
medx2=0.1374
varx2=0.0024
sdx2=varx2^0.5
cov_x1x2=-cor_x1x2*(sdx1*sdx2)

sda=(varx1-cov_x1x2)^0.5
sdb=(varx2-cov_x1x2)^0.5
a=rnorm(n,medx1,sda);var(a)
b=rnorm(n,medx2,sdb);var(b)
c=rnorm(n,0,cov_x1x2^0.5);var(c)

x1=a+c
x2=b-c
var(x1)
var(x2)
cor(x1,x2)


Abraços,

Luiz Roberto.


Luiz Roberto Martins Pinto
Prof. Pleno/DCET/UESC
Laboratório de Estatística Computacional
Universidade Estadual de Santa Cruz
Ilhéus-Bahia-Brasil

luizroberto.u...@gmail.com
skype: lrmpinto
http://lattes.cnpq.br/2732314327604831

"*The s**cience exists because there are patterns. The patterns exist
because God created them*.
* The statistic exists to research the patterns that God created.*"



Em 18 de outubro de 2016 14:52, Éder Comunello via R-br <
r-br@listas.c3sl.ufpr.br> escreveu:

> Adriele, boa tarde!
>
> Tentei adaptar seu código R com base no que entendi do script em SAS. Veja
> se pode ser útil...
>
> ### <code r>
> library(MASS)
> {
>   set.seed(1357)
>   n  <- 1000  # tamanho da amostra gerada (N)
>   p  <- 2     # numero de variáveis a serem geradas (K)
>   ME     <- rep(1, p)
>   rho    <- -0.2 # correlacao negativa - 0.2
>   sigma2 <- 1
>   sigma  <- sigma2 * ((1-rho)*diag(p)+rho*matrix(1, p, p))
>   SI  <- mvrnorm(n, ME, sigma)
> }
> apply(SI, 2, mean) # ~ 1
> # [1] 1.037774 1.017846
>
> cor(SI)            # ~ sigma
> # [,1]       [,2]
> # [1,]  1.0000000 -0.1973338
> # [2,] -0.1973338  1.0000000
>
> mu <- c(0.0067, 0.1374)
> mean(0.0017*SI[,1]) # ~ 0.0017
> mean(0.0024*SI[,2]) # ~ 0.0024
>
> RES <- data.frame(X1=0.0067+0.0017*SI[,1], X2=0.1374+0.0024*SI[,2])
>
> round(cbind(RES=colMeans(RES), INI=mu),4)
> # RES    INI
> # X1 0.0085 0.0067
> # X2 0.1398 0.1374
>
> cor(RES)
> # X1         X2
> # X1  1.0000000 -0.1973338
> # X2 -0.1973338  1.0000000
> ### </code>
>
>
> ================================================
> Éder Comunello
> Researcher at Brazilian Agricultural Research Corporation (Embrapa)
> DSc in Agricultural Systems Engineering (USP/Esalq)
> MSc in Environ. Sciences (UEM), Agronomist (UEM)
> ---
> Embrapa Agropecuária Oeste, Dourados, MS, Brazil |<O>|
> ================================================
> GEO, -22.2752, -54.8182, 408m
> UTC-04:00 / DST: UTC-03:00
>
>
>
>
> Em 17 de outubro de 2016 16:05, Adriele Giaretta Biase via R-br <
> r-br@listas.c3sl.ufpr.br> escreveu:
>
>> Olá pessoal,
>>
>>
>> tenho uma dúvida com relação à geração de variáveis aleatórias normais
>> multivariadas no R.
>>
>> Eu gostaria de gerar duas variáveis aleatórias (x1, x2) usando
>> distribuição normal multivariada com estrutura de correlação entre elas. A
>> Estrutura da matriz de correlação foi estimada antes com base num banco de
>> dados reais (com correlação fraca de - 0.2). Porém, essas variáveis não
>> podem ser negativas, Ex. x1 tem média e variância, respectivamente, iguais
>> a 0.0067 e 0.0017; x2 tem média e variância, respectivamente, iguais a
>> 0.1374 e 0.0024. Eu gero as variáveis da seguinte forma:
>>
>>
>>
>> *Programa usado no R:*
>>
>> # Simulando os parâmetros com estrutura de correlação entre eles
>>
>> n <- 1000  # tamanho da amostra gerada
>>
>> p <- 2  # numero de variáveis a serem geradas
>>
>>
>>
>> library(MASS)
>>
>> # construíndo a matriz de correlação para usar nas simulações, baseadas
>> na característica da amostra coletada
>>
>> mu<-rep(0, times = p)
>>
>> rho <-   - 0.2             # correlacao negativa - 0.2
>>
>> sigma2 <- 1
>>
>> Sigma <- sigma2 * ((1-rho)*diag(p)+rho*matrix(1, p, p))
>>
>> X1 <- 0.0096  # Media da variavel x1
>>
>> X2 <- 0.1203 # Media da variavel x2
>>
>> media <- c(X1,X2)
>>
>> y  <-  mvrnorm(n, media, Sigma)
>>
>> cor(y)
>>
>> apply(y, 2, mean)
>>
>> y
>>
>>
>>    [1,]  0.309910452  1.0642521083
>>
>>    [2,] -0.251583312  1.8909032058
>>
>>    [3,]  1.330362012 -1.1239501814
>>
>>    [4,] -0.793399464 -1.5433056284
>>
>>    [5,]  2.165144843 -0.2645184534
>>
>>    [6,]  0.532777085  1.0910864562
>>
>>    [7,] -1.612135390  2.0489354648
>>
>>    [8,] -0.430529913  1.0312062602
>>
>> ...
>>
>>
>>
>> Quando gero as simulações para x1 e x2 usando a função *mvrnorm*, o
>> resultado me retorna alguns valores negativos para as variáveis, isso não
>> poderia ocorrer. Teria alguma outra função em que eu possa fornecer a
>> variância de cada uma dessas variáveis, além da estrutura de correlação?
>> Como poderia contornar essa situação, usando o R?
>>
>>
>>
>> *Obs:* No SAS, a seguinte programação funciona perfeitamente:
>>
>>
>> *proc* *iml*;
>>
>> wrksize=*100000*;
>>
>> K=*2*;
>>
>> N=*1000*; /* tamanho da amostra */
>>
>> M={*0* *0*};
>>
>> S={ *1*       -*0.5283*,
>>
>>    -*0.5283*       *1* };
>>
>> X=shape(*0*,K,N);
>>
>> ME=*0*; SI=*1*;
>>
>> DO I=*1* TO K;
>>
>>   DO J=*1* TO N;
>>
>>     if I>*1*
>>
>>       then
>>
>>         do;
>>
>>          ME=M[I]+(S[*1*:I-*1*,I])`*(inv(S[*1*:I-*1*,*1*:I-*1*])*(X[*1*:I-
>> *1*,J]-M[*1*:I-*1*]));
>>
>>          SI=S[I,I]-(S[*1*:I-*1*,I])`*(inv(S[*1*:I-*1*,*1*:I-*1*])*(S[*1*
>> :I-*1*,I]));
>>
>>         end;
>>
>>     X[I,J]=ME+NORMAL(*0*)*SQRT(SI);
>>
>>
>>
>>   END;
>>
>>
>>
>> END;
>>
>> Z=t(X);
>>
>>
>>
>> varnames='X1':'X2';
>>
>> create NOVO from Z[colname=varnames];
>>
>> append from Z;
>>
>> *quit*;
>>
>>
>>
>> *data* MEXT; set NOVO;
>>
>> options ps=*66* ls=*75*;
>>
>>    X1=*0.0067*+*0.0017**X1; /* normal */
>>
>>    X2=*0.1374*+*0.0024**X2; /* normal */
>>
>>
>>
>> *proc* *corr* data=MEXT;
>>
>> var X1 X2;
>>
>> *run*;
>>
>> --
>> Adriele Giaretta Biase.
>> Mestre em  Estatística e Experimentação Agropecuária - UFLA.
>> Doutora em Estatística e Experimentação Agronômica - ESALQ/ USP
>> Contato: (19) 98861-0619.
>>
>> _______________________________________________
>> R-br mailing list
>> R-br@listas.c3sl.ufpr.br
>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
>> código mínimo reproduzível.
>>
>
>
> _______________________________________________
> R-br mailing list
> R-br@listas.c3sl.ufpr.br
> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
> código mínimo reproduzível.
>
_______________________________________________
R-br mailing list
R-br@listas.c3sl.ufpr.br
https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forne�a c�digo 
m�nimo reproduz�vel.

Responder a