模擬分佈圖

Question

我正在研究具有不同參數（樣本大小和方差）的正態分佈和伽馬分佈的魯棒性。

我也得到了模擬結果。他們是三張桌子。 但現在我不得不嘗試繪製模擬分佈圖，以使人們更好地理解結果。

我在R中還是個新手。我是否需要在分佈圖中包含所有三個表的結果？

######################################################################## 
#For gamma distribution with equal skewness 1.5 

# rm(list=ls()) # clean the workspace 
nSims<-10000  #set the number of simulations 
alpha<-0.05  #set the significance level 

# to ensure the reproduction of the result 
# here we declare the random seed generator 
set.seed(1) 

#create vector to combine all std deviations 
sds<-matrix(c(4,4,6,4,8,4,10,4,12,4),nrow=2) 

sd1<-c(4,6,8,10,12) 
sd2<-c(4,4,4,4,4) 

## Put the samples sizes into matrix then use a loop for sample sizes 
sample_sizes<-matrix(c(10,10,10,25,25,25,25,50,25,100,50,25,50,100,100,25,100,100), 
nrow=2) 

#shape parameter for gamma distribution for equal skewness 
#forty five cases for each skewness!!!! 
sp1<-matrix(rep(c(16/9),each=45),ncol=1) 

scp <- c(1,1.5,2,2.5,3) 

##(use expand.grid)to create a data frame 
ss_scp<- expand.grid(sample_sizes[2,],scp) 

#create a matrix combining the forty five cases of combination of sample sizes,shape and scale parameter 
all <- cbind(rep(sample_sizes[1,], 5),ss_scp[,1],sp1,ss_scp[,2]) 

# name the column samples 1 and 2 and standard deviation 
colnames(all) <- c("m","n","sp","scp") 

#set empty vector of length no.of simulation(10000) to store p-value 
equal<-unequal<-mann<-c(rep(0,nrow(all))) 

#set nrow =nsims because wan storing every p-value simulated 
#for gamma distribution with equal skewness 
matrix_t <-matrix(0,nrow=nSims,ncol=5) 
matrix_u<-matrix(0,nrow=nSims,ncol=5) 
matrix_mann <-matrix(0,nrow=nSims,ncol=5) 

##for the samples sizes into matrix then use a loop for sample sizes 
# this loop steps through the all_combine matrix 
    for(ss in 1:nrow(all)) 
    { 
    #generate samples from the first column and second column 
    m<-all[ss,1] 
    n<-all[ss,2] 

     for (sim in 1:nSims) 
     { 
     #generate 2 random samples from gamma distribution with equal skewness 
     gamma1<-rgamma(m,all[ss,3],scale=all[ss,4]) 
     gamma2<-rgamma(n,all[ss,3],scale=1) 

     gamma1<-gamma1-all[ss,3]*all[ss,4] 
     gamma2<-gamma2-all[ss,3] 

     #extract p-value out and store every p-value into matrix 
     p<-t.test(gamma1,gamma2,var.equal=TRUE)$p.value 
     q<-t.test(gamma1,gamma2,var.equal=FALSE)$p.value 
     r<-wilcox.test(gamma1,gamma2)$p.value 

     matrix_t[sim,1]<- p 
     matrix_u[sim,1]<- q 
     matrix_mann[sim,1] <- r 
    } 
     ##store the result 
     equal[ss]<- sum(matrix_t[,1]<alpha) 
     unequal[ss]<-sum(matrix_u[,1]<alpha) 
     mann[ss]<- sum(matrix_mann[,1]<alpha) 
    } 

g1_equal<-cbind(all, equal, unequal, mann) 
print("g1_equal_skewness1.5)") 
print(g1_equal) 

    #samples sizes (10,10),(10,25).. 
    #standard deviation ratio (1,1.5,2,2.5,3) 
        Gamma(equal skewness)   Gamma(unequal skewness) 
        1.5 2.0 2.5 3.0 3.5  (1.5,1) 2,1.5 2.5,2 3,2.5 3.5,3 
10,10 
      Normal 
    1.0  506  382 379 343 270  246  422 426 383 303 247 
    1.5  472  493 463 507 537  571  531 518 548 528 532 
    2.0 516  597 679 736 829  935  597 680 760 836 951 
    2.5 498  627 747 905 1028 1215  687 825 944 1011 1197 
    3.0 493  678 864 1010 1190 1379  705 831 1015 1170 1436 


10,25 

    1.0 511  568 557 633 647 630  603 599 604  652 654 
    1.5 501  692 840 977 1012 1173  675 756 940  1068 1130 
    2.0 438  713 951 1049 1264 1470  773 869 1055 1259 1401 
    2.5 506  810 939 1101 1300 1594  761 960 1155 1339 1512 
    3.0 524  787 933 1176 1378 1599  772 967 1201 1339 1612 



25,25 

    1.0 479  463 451 447 417 414  513 429 439 469 392 
    1.5 493  534 556 504 568 587  537 517 528 539 555 
    2.0 510  543 599 676 663 773  538 607 677 712 725 
    2.5 487  591 662 731 807 908  581 643 733 769 893 
    3.0 488  614 668 761 811 1002 582 694 728 900 946 


25,50 

    1.0 519  585 487 569 559 579  521 572 568 581 583 
    1.5 510  532 651 695 725 836  625 647 729 737 802 
    2.0 501  586 660 758 846 888  618 653 794 876 957 
    2.5 466  635 687 823 937 996  612 702 782 909 1025 
    3.0 492  603 719 824 970 1045 640 704 826 945 1073 


25,100 

    1.0 486  559 589 670 726 778 552 614 666 752 750 
    1.5 494  621 700 787 903 955 602 703 774 842 1008 
    2.0 516  617 707 817 969 1073 613 755 774 932 1091 
    2.5 470  598 731 873 969 1118 624 752 849 970 1094 
    3.0 493  710 718 824 1021 1167 645 746 887 988 1149 

50,25 

    1.0 495  507 511 552 550 534  491 527 496 554 534 
    1.5 535  472 470 489 470 413  458 503 460 456 410 
    2.0 499  507 478 488 468 465  495 490 542 528 489 
    2.5 486  500 532 517 559 629  509 493 526 569 601 
    3.0 490  586 536 561 654 644  544 567 563 614 665 

50,100 1.0 518  515 530 531 514 569  516 494 517 548 578 
     1.5 528  503 542 597 596 656  554 565 612 606 708 
     2.0 453  525 588 640 727 775  520 625 628 727 772 
     2.5 500  586 660 669 733 837  552 622 660 695 802 
     3.0 494  557 640 680 747 847  582 634 686 776 834 

100,25             
     1.0 489  553 607 641 712 777  557 560 653 677 751 
     1.5 516  497 553 532 619 595  496 548 512 549 553 
     2.0 500  492 483 518 472 468  536 521 497 463 463 

     2.5 493  498 473 446 488 461  483 463 476 452 472 
     3.0 482  490 516 481 488 500  563 477 496 492 537 




100,100 
     1.0 472  508 492 483 517 487  517 521 476 505 485 
     1.5 507  498 496 511 518 546  520 520 498 547 531 
     2.0 465  478 540 542 584 599  496 504 585 558 589 
     2.5 508  486 566 551 614 602  520 539 583 601 642 
     3.0 494  497 575 545 614 651  561 557 590 615 624 

Answer 1

這聽起來像你想知道如何着手繪製你的分佈。您可以創建三個獨立的圖，也可以將每個分佈覆蓋在一個圖上。以下是如何在ggplot2中做到這一點。

假設您的數據框爲df，列dist1,dist2和dist3。

install.packages('ggplot2') 
library(ggplot2) 

ggplot(df, aes(x = dist1)) + 
    geom_density(color = 'red') + 
    geom_density(aes(x = dist2), color = 'green') + 
    geom_density(aes(x = dist3), color = 'blue') 

這應該給你一個三行密度圖，每個分佈一個。如果你想創建三個獨立的地塊，只需爲每個分配創建一個新的ggplot。

plot1 <- ggplot(df, aes(x = dist1) + geom_density() 
plot2 <- ggplot(df, aes(x = dist2) + geom_density() 

...等等。這有幫助嗎？