2016-03-24 80 views
0

我的代碼存儲在矩陣的值:使用循環

set.seed(100) 
n <- 10 
d <- rep(NA, n) 
d[1] <- 0 
y <- runif(n) 
a <- 10 

for (i in 2:length(y)) { 
    d[i] <- d[i-1] + y[i-1] 
} #This creates an interval with endpoints at my y uniform RVs 

store.x <- NULL 
for (j in 1:a) { 
    x <- runif(1, min =0, max =sum(y)) 
    for (i in 1:length(y)) { 
     if (x <= d[i+1] && x > d[i]) { 
      store.x[j] <- i 
      break 
     } 
    } 
} #This tells you which interval my x uniform RV is in 

現在不是store.x是,告訴我什麼時間間隔X落入一個載體,我希望它被存儲在相應的行和列的值爲1.因此,對於我的第一個x,因爲它落入區間7,我的矩陣將全部爲零,除了第一行,第七列和第七行中的一個,第一個柱。

關於如何做到這一點的任何想法將不勝感激! 感謝

+0

爲什麼要「除了第一行第七列和第七行第一列中的一個都是零」?是不是「除了第一行第七列中的所有零」都已經足夠了? – adaien

+0

@adiana就這樣對稱 – Killian

+0

所以a應該總是與長度(y)相同? – adaien

回答

1

可以使用矢量調用runif()cumsum()更有效地實現這個算法,並findInterval()

set.seed(1L); ## seed the PRNG for reproducible results 
n <- 10L; ## length of x, y, and result matrix dimensions 
y <- runif(n); ## produce random interval lengths 
yb <- cumsum(c(0,y)); ## calculate interval boundaries 
x <- runif(n,0,yb[length(yb)]); ## generate random x values 
i <- findInterval(x,yb); ## find which intervals contain the x values 
m <- matrix(0,n,n); ## init the result matrix 
m[matrix(c(seq_len(n),i,i,seq_len(n)),ncol=2L)] <- 1; ## symmetric 1 assignment 
m; ## print result 
##  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] 
## [1,] 0 0 1 0 0 0 0 0 0  0 
## [2,] 0 0 1 0 0 0 0 0 0  0 
## [3,] 1 1 0 0 0 0 1 0 0  0 
## [4,] 0 0 0 1 0 0 0 0 1  0 
## [5,] 0 0 0 0 0 0 0 1 0  0 
## [6,] 0 0 0 0 0 1 0 0 0  0 
## [7,] 0 0 1 0 0 0 1 0 0  0 
## [8,] 0 0 0 0 1 0 0 0 0  1 
## [9,] 0 0 0 1 0 0 0 0 0  0 
## [10,] 0 0 0 0 0 0 0 1 0  0 

從技術上說,你寫的間隔測試x <= d[i+1] && x > d[i]的方式意味着你要左邊界爲,打開,右邊界爲,每個區間關閉。看起來findInterval()最近纔在其邏輯上增加了對此變體的支持,特別是在r69814開發快照中(最終將變爲R-3.3.0)。所以你可能還沒有訪問它,但是當你這樣做時,你可以通過left.open=T來獲得該行爲。


如果我們有xa元件允許a != n那麼我們必須決定如何各元素的索引中x映射到結果矩陣的索引。上述解決方案假定a == n並且兩個索引域之間存在直接映射。

如果我們考慮的x元件以對應於導致矩陣索引1:n以循環方式,那麼我們可以定義xi到是循環索引,取yii,即其中x元件降落y間隔指數。然後,我們可以彙總所有(xi,yi)命中以產生每個單元格的計數。

如果你仍然需要對稱性,那麼我們可以進一步累積每個(yi,xi)的計數,從而對每個擊中的每個擊打進行重複計數,每個擊中兩個對稱細胞。

set.seed(1L); ## seed the PRNG for reproducible results 
n <- 10L; ## length of y and result matrix dimensions 
a <- 100L; ## length of x 
y <- runif(n); ## produce random interval lengths 
yb <- cumsum(c(0,y)); ## calculate interval boundaries 
x <- runif(a,0,yb[length(yb)]); ## generate random x values 
i <- findInterval(x,yb); ## find which intervals contain the x values 
m <- matrix(0,n,n); ## init the result matrix 
hit <- as.matrix(aggregate(n~xi+yi,cbind(xi=seq_len(n),yi=i,n=1),sum)); ## aggregate hits 
m[hit[,c('xi','yi')]] <- m[hit[,c('xi','yi')]]+hit[,'n']; ## add n to xi,yi 
m[hit[,c('yi','xi')]] <- m[hit[,c('yi','xi')]]+hit[,'n']; ## add n to yi,xi 
m; ## print result 
##  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] 
## [1,] 0 0 1 2 0 3 3 1 2  0 
## [2,] 0 2 2 2 1 1 3 3 1  0 
## [3,] 1 2 0 4 1 5 4 2 0  1 
## [4,] 2 2 4 10 1 2 1 1 3  2 
## [5,] 0 1 1 1 0 3 2 6 0  2 
## [6,] 3 1 5 2 3 2 3 5 1  0 
## [7,] 3 3 4 1 2 3 4 2 3  3 
## [8,] 1 3 2 1 6 5 2 2 3  2 
## [9,] 2 1 0 3 0 1 3 3 2  2 
## [10,] 0 0 1 2 2 0 3 2 2  0 
+0

這太棒了,謝謝!如果我想爲作業添加1而不是僅僅調用1,該怎麼辦?所以如果我有幾百個,我可以看到有多少行/列 – Killian

+0

@Killian請參閱編輯。 – bgoldst

1
store.x <- matrix(0,nrow=length(y),ncol=length(y)) 
for(j in 1:length(y)) { 
x <- runif(1, min= 0, max =sum(y)) 
for(i in 1:length(y)) { 
    if(x <= d[i+1] && x > d[i]) { 
     store.x[j,i] <- 1 
     store.x[i,j] <- 1 
     break 
    }}}