Bivariate Normal Distribution Plots

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|

Bivariate Normal Distribution Plots

reichmaj
R-Help

I am attempting to create a series of bivariate normal distributions.  So using the mvtnorm library I have created the following code ...

# Standard deviations and correlation
sig_x <- 1
sig_y <- 1
rho_xy <- 0.0

# Covariance between X and Y
sig_xy <- rho_xy * sig_x *sig_y

# Covariance matrix
Sigma_xy <- matrix(c(sig_x ^ 2, sig_xy, sig_xy, sig_y ^ 2), nrow = 2, ncol = 2)

# Load the mvtnorm package
library("mvtnorm")

# Means
mu_x <- 0
mu_y <- 0

# Simulate 1000 observations
set.seed(12345)  # for reproducibility
xy_vals <- rmvnorm(1000, mean = c(mu_x, mu_y), sigma = Sigma_xy)

# Have a look at the first observations
head(xy_vals)

# Create scatterplot
plot(xy_vals[, 1], xy_vals[, 2], pch = 16, cex = 2, col = "blue",
     main = "Bivariate normal: rho = 0.0", xlab = "x", ylab = "y")

# Add lines
abline(h = mu_y, v = mu_x)

Problem is this results in sigma(x) = sigma(y), rho=0 and I need or what 2sigma(x)=sigma(y), rho=0 or 2sigma(y)=sigma(x), rho=0 to elongate the distribution.  What I have created creates a circle.  Can I do that within the mvtnorm package?

Jeff Reichman

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Bivariate Normal Distribution Plots

S Ellison-2


> -----Original Message-----
> From: R-help [mailto:[hidden email]] On Behalf Of JEFFERY REICHMAN


> # Standard deviations and correlation
> sig_x <- 1
> sig_y <- 1
> rho_xy <- 0.0
>
> # Covariance between X and Y
> sig_xy <- rho_xy * sig_x *sig_y
>
> # Covariance matrix
> Sigma_xy <- matrix(c(sig_x ^ 2, sig_xy, sig_xy, sig_y ^ 2), nrow = 2, ncol = 2)
> ...
> Problem is this results in sigma(x) = sigma(y), rho=0 ...
> What I have created creates a circle.  
Er, yes ... that is what you asked for.

Have you tried rho_xy=0.7 or similar in the code above?


*******************************************************************
This email and any attachments are confidential. Any use...{{dropped:8}}

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Bivariate Normal Distribution Plots

Richard M. Heiberger
In reply to this post by reichmaj
Please look at my book
Statistical Analysis and Data Display
https://www.springer.com/us/book/9781493921218

Figures 3.8, 3.9, 3.10

The code for these figures is available in the HH package

install.packages("HH")
library(HH)
HHscriptnames(3) ## this gives the filename on your computer containing the code
## open the file in your preferred editor and run chunks 15 and 16


On Thu, Apr 12, 2018 at 10:59 AM, JEFFERY REICHMAN
<[hidden email]> wrote:

> R-Help
>
> I am attempting to create a series of bivariate normal distributions.  So using the mvtnorm library I have created the following code ...
>
> # Standard deviations and correlation
> sig_x <- 1
> sig_y <- 1
> rho_xy <- 0.0
>
> # Covariance between X and Y
> sig_xy <- rho_xy * sig_x *sig_y
>
> # Covariance matrix
> Sigma_xy <- matrix(c(sig_x ^ 2, sig_xy, sig_xy, sig_y ^ 2), nrow = 2, ncol = 2)
>
> # Load the mvtnorm package
> library("mvtnorm")
>
> # Means
> mu_x <- 0
> mu_y <- 0
>
> # Simulate 1000 observations
> set.seed(12345)  # for reproducibility
> xy_vals <- rmvnorm(1000, mean = c(mu_x, mu_y), sigma = Sigma_xy)
>
> # Have a look at the first observations
> head(xy_vals)
>
> # Create scatterplot
> plot(xy_vals[, 1], xy_vals[, 2], pch = 16, cex = 2, col = "blue",
>      main = "Bivariate normal: rho = 0.0", xlab = "x", ylab = "y")
>
> # Add lines
> abline(h = mu_y, v = mu_x)
>
> Problem is this results in sigma(x) = sigma(y), rho=0 and I need or what 2sigma(x)=sigma(y), rho=0 or 2sigma(y)=sigma(x), rho=0 to elongate the distribution.  What I have created creates a circle.  Can I do that within the mvtnorm package?
>
> Jeff Reichman
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Bivariate Normal Distribution Plots

R help mailing list-2
In reply to this post by reichmaj
Try this code:

# Standard deviations and correlation
sig_x <- 1
sig_y <- 2
rho_xy <- 0.7

# Covariance between X and Y
sig_xy <- rho_xy * sig_x *sig_y

# Covariance matrix
Sigma_xy <- matrix(c(sig_x ^ 2, sig_xy, sig_xy, sig_y ^ 2), nrow = 2,
ncol = 2)

# Load the mvtnorm package
library("mvtnorm")

# Means
mu_x <- 0
mu_y <- 0

# Simulate 1000 observations
set.seed(12345)  # for reproducibility
xy_vals <- rmvnorm(1000, mean = c(mu_x, mu_y), sigma = Sigma_xy)

# Have a look at the first observations
head(xy_vals)

# Create scatterplot
# plot(xy_vals[, 1], xy_vals[, 2], pch = 16, cex = 2, col = "blue",
#      main = "Bivariate normal: rho = 0.0", xlab = "x", ylab = "y")

library(graphics)

x <- xy_vals[, 1]
y <- xy_vals[, 2]

par(mar=c(4, 4, 2, 6)+0.4)

smoothScatter(x, y, asp=1,
               main = paste("Bivariate normal: rho = ", rho_xy),
               xlab = "x", ylab = "y")


# Add lines
abline(h = mu_y, v = mu_x)

library(fields)

n <- matrix(0, ncol=128, nrow=128)

xrange <- range(x)
yrange <- range(y)

for (i in 1:length(x)) {
   posx <- 1+floor(127*(x[i]-xrange[1])/(xrange[2]-xrange[1]))
   posy <- 1+floor(127*(y[i]-yrange[1])/(yrange[2]-yrange[1]))
   n[posx, posy] <- n[posx, posy]+1
}

image.plot( legend.only=TRUE,
             zlim= c(0, max(n)), nlevel=128,
             col=colorRampPalette(c("white", blues9))(128))


Hope it helps,

Marc

Le 12/04/2018 à 16:59, JEFFERY REICHMAN a écrit :

> R-Help
>
> I am attempting to create a series of bivariate normal distributions.  So using the mvtnorm library I have created the following code ...
>
> # Standard deviations and correlation
> sig_x <- 1
> sig_y <- 1
> rho_xy <- 0.0
>
> # Covariance between X and Y
> sig_xy <- rho_xy * sig_x *sig_y
>
> # Covariance matrix
> Sigma_xy <- matrix(c(sig_x ^ 2, sig_xy, sig_xy, sig_y ^ 2), nrow = 2, ncol = 2)
>
> # Load the mvtnorm package
> library("mvtnorm")
>
> # Means
> mu_x <- 0
> mu_y <- 0
>
> # Simulate 1000 observations
> set.seed(12345)  # for reproducibility
> xy_vals <- rmvnorm(1000, mean = c(mu_x, mu_y), sigma = Sigma_xy)
>
> # Have a look at the first observations
> head(xy_vals)
>
> # Create scatterplot
> plot(xy_vals[, 1], xy_vals[, 2], pch = 16, cex = 2, col = "blue",
>       main = "Bivariate normal: rho = 0.0", xlab = "x", ylab = "y")
>
> # Add lines
> abline(h = mu_y, v = mu_x)
>
> Problem is this results in sigma(x) = sigma(y), rho=0 and I need or what 2sigma(x)=sigma(y), rho=0 or 2sigma(y)=sigma(x), rho=0 to elongate the distribution.  What I have created creates a circle.  Can I do that within the mvtnorm package?
>
> Jeff Reichman
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.