# Efficient code

**URL:** <https://discourse.julialang.org/t/efficient-code/5508>\
**Category:** General Usage\
**Tags:** performance\
**Created:** [August 22, 2017, 4:22pm UTC](https://discourse.julialang.org/t/efficient-code/5508 "2017-08-22T16:22:34Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![bdeonovic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bdeonovic/32/3928_2.png) [@bdeonovic](https://discourse.julialang.org/u/bdeonovic)\
**Post date:** [August 22, 2017, 4:22pm UTC](https://discourse.julialang.org/t/efficient-code/5508/1 "2017-08-22T16:22:34Z")

</div>

Setup for MWE

```julia
using Distributions
qidx = collect(1:10)
x = rand([0,1], 10)
theta = randn()
gamma=randn(2)
d = randn(10)
A = [.3 .7; .2 .8]
p=[.8, .2]

```

I’d like to run the following code, but I’m worried whether I’ve written this efficiently. For my example `m` will stay as `2` but `n` will be much larger.

This is forward-backward algo for HMMs in matrix notation but with a time varying emission matrix.

```julia
 m = length(p)                                                      
 n = length(x)                                                      
 G = zeros(n, m)                                                    
 l = zeros(m)                                                       
 c = zeros(n)                                                       
                                                                    
 B = zeros(m, m, n)                                                 
  for j in 1:n                                      
   s = cdf(Logistic(theta - gamma[2]), d[qidx[j]]) 
   g = cdf(Logistic(d[qidx[j]] - gamma[1]), theta) 
   if x[j] == 1                                    
     B[:, :, j] = diagm([1.0 - s, g])              
   else                                            
     B[:, :, j] = diagm([s, 1.0 - g])              
   end                                             
 end                                               
                                                   
 G[1,:] = B[:,:,1]'p                              
                                                   
 for j in 2:n                                      
   G[j,:] = G[j-1,:]' * A * B[:,:,j]             
 end                                               

```

In particular I am interested in whether the way I have allocated the `B` matrix is good (the last dimension is what is looped over), would it be better if it was the first dimension? I always get confused with which dimension should be looped over first, second, etc.

* * *

EDIT:  
I’ve changed `B` so that the first dimension is used in the innermost loop.

```julia
 B = zeros(n, m, m)                                     
                                                        
 for j in 1:n                                           
   s = cdf(Logistic(theta - gamma[2]), d[qidx[j]])      
   g = cdf(Logistic(d[qidx[j]] - gamma[1]), theta)      
   if x[j] == 1                                         
     B[j,:, :] = diagm([1.0 - s, g])                   
   else                                                 
     B[j,:,:] = diagm([s, 1.0 - g])                    
   end                                                  
 end                                                    
                                                        
 G[1,:] = B[1,:,:]'p                                  
                                                        
 for j in 2:n                                           
   G[j,:] = G[j-1,:]' * A * B[j,:,:]                 
 end                                                                                 

```

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [August 22, 2017, 4:31pm UTC](https://discourse.julialang.org/t/efficient-code/5508/2 "2017-08-22T16:31:22Z")

</div>

> [@bdeonovic](#):
>
> In particular I am interested in whether the way I have allocated the B matrix is good (the last dimension is what is looped over), would it be better if it was the first dimension? I always get confused with which dimension should be looped over first, second, etc.

The first dimension should be used for innermost loops with Julia’s column-major storage.