Kolumnowa normalizacja macierzy w R [zamknięty]

26

Chciałbym przeprowadzić kolumnową normalizację macierzy w R. Biorąc pod uwagę macierz m, chcę znormalizować każdą kolumnę dzieląc każdy element przez sumę kolumny. Jeden (hackish) sposób, aby to zrobić, jest następujący:

m / t(replicate(nrow(m), colSums(m)))

Czy istnieje bardziej zwięzły / elegancki / skuteczny sposób na osiągnięcie tego samego zadania?

mavam
źródło

Odpowiedzi:

41

Po to są zamiatanie i skala.

sweep(m, 2, colSums(m), FUN="/")
scale(m, center=FALSE, scale=colSums(m))

Alternatywnie możesz użyć recyklingu, ale musisz go transponować dwukrotnie.

t(t(m)/colSums(m))

Lub możesz skonstruować pełną macierz, którą chcesz podzielić, tak jak w swoim pytaniu. Oto inny sposób, w jaki możesz to zrobić.

m/colSums(m)[col(m)]

Zwróć też uwagę na dodatek karakala z komentarzy:

m %*% diag(1/colSums(m))
Aaron - Przywróć Monikę
źródło
8
Jeszcze jedno:m %*% diag(1/colSums(m))
karakal
Nigdy wcześniej nie słyszałem o funkcji zamiatania, dzięki!
Matteo De Felice,
10

Innym jest prop.table(m, 2), lub po prostu propr(m), wewnętrzny sweep.

Może być interesujące porównanie wydajności tych równoważnych rozwiązań, więc zrobiłem mały test porównawczy (używając microbenchmarkpakietu).

Oto matryca wejściowa, mktórej użyłem:

          [,1]         [,2]         [,3]         [,4]         [,5]
A 1.831564e-02 4.978707e-02 1.353353e-01 3.678794e-01 3.678794e-01
B 3.678794e-01 1.353353e-01 4.978707e-02 1.831564e-02 6.737947e-03
C 4.539993e-05 2.061154e-09 9.357623e-14 4.248354e-18 5.242886e-22
D 1.831564e-02 4.978707e-02 1.353353e-01 3.678794e-01 3.678794e-01
E 3.678794e-01 1.353353e-01 4.978707e-02 1.831564e-02 6.737947e-03
F 4.539993e-05 2.061154e-09 9.357623e-14 4.248354e-18 5.242886e-22
G 1.831564e-02 4.978707e-02 1.353353e-01 3.678794e-01 3.678794e-01
H 3.678794e-01 1.353353e-01 4.978707e-02 1.831564e-02 6.737947e-03
I 4.539993e-05 2.061154e-09 9.357623e-14 4.248354e-18 5.242886e-22

Oto konfiguracja testu porównawczego:

microbenchmark(
prop = prop.table(m, 2),
scale = scale(m, center=FALSE, scale=colSums(m)),
sweep = sweep(m, 2, colSums(m), FUN="/"),
t_t_colsums = t(t(m)/colSums(m)),
m_colsums_col = m/colSums(m)[col(m)],
m_mult_diag = m %*% diag(1/colSums(m)),
times = 1500L)

Oto wyniki testu porównawczego:

Unit: microseconds
           expr     min       lq   median       uq      max
1 m_colsums_col  29.089  32.9565  35.9870  37.5215 1547.972
2   m_mult_diag  43.278  47.6115  51.7075  53.8945  110.560
3          prop 207.070 214.3010 216.6800 219.9680 2091.913
4         scale 133.659 142.6325 145.3100 147.9195 1730.640
5         sweep 113.969 119.6315 121.3725 123.6570 1663.356
6   t_t_colsums  56.976  65.3580  67.8895  69.5130 1640.660

Dla kompletności jest to wynik:

          [,1]         [,2]         [,3]         [,4]         [,5]
A 1.580677e-02 8.964714e-02 2.436862e-01 3.175247e-01 3.273379e-01
B 3.174874e-01 2.436862e-01 8.964714e-02 1.580862e-02 5.995403e-03
C 3.918106e-05 3.711336e-09 1.684944e-13 3.666847e-18 4.665103e-22
D 1.580677e-02 8.964714e-02 2.436862e-01 3.175247e-01 3.273379e-01
E 3.174874e-01 2.436862e-01 8.964714e-02 1.580862e-02 5.995403e-03
F 3.918106e-05 3.711336e-09 1.684944e-13 3.666847e-18 4.665103e-22
G 1.580677e-02 8.964714e-02 2.436862e-01 3.175247e-01 3.273379e-01
H 3.174874e-01 2.436862e-01 8.964714e-02 1.580862e-02 5.995403e-03
I 3.918106e-05 3.711336e-09 1.684944e-13 3.666847e-18 4.665103e-22

Bez wątpienia wygrywa mała matryca !m / colSums(m)[col(m)]


Ale dla dużych matryc? W kolejnym przykładzie użyłem matrycy 1000 x 1000.

set.seed(42)
m <- matrix(sample(1:10, 1e6, TRUE), 1e3)
...
Unit: milliseconds
           expr      min       lq   median        uq       max
1 m_colsums_col 55.26442 58.94281 64.41691 102.69683 119.08685
2   m_mult_diag 34.67692 41.68494 80.05480  89.48099  99.72062
3          prop 87.95552 94.13143 99.17044 136.03669 160.51586
4         scale 52.84534 55.07107 60.57154  99.87761 156.16622
5         sweep 52.79542 55.93877 61.55066  99.67766 119.05134
6   t_t_colsums 63.09783 65.53783 68.93731 110.03691 127.89792

Dla dużych matryc m / colSums(m)[col(m)] działa dobrze (4. pozycja), ale nie wygrywa .

Dla dużych matryc m %*% diag(1/colSums(m)) wygrywa !

leodido
źródło
1
jaki pakiet jest proprz?
Glen_b
5
apply(m,2,norm<-function(x){return (x/sum(x)}) ?
Sowmya Iyer
źródło
4
Witamy na stronie @Sowmyalyer. Czy miałbyś coś przeciwko dodaniu tekstu w celu pełniejszego przedstawienia i wyjaśnienia swojej odpowiedzi?
gung - Przywróć Monikę