Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Comments on stdlib_stats_corr.f90 #955

Open
Beliavsky opened this issue Mar 18, 2025 · 2 comments
Open

Comments on stdlib_stats_corr.f90 #955

Beliavsky opened this issue Mar 18, 2025 · 2 comments

Comments

@Beliavsky
Copy link

Beliavsky commented Mar 18, 2025 •
edited
Loading

I think there are a number of flaws in corr_mask_2_rsp_rsp from https://github.com/fortran-lang/stdlib/blob/stdlib-fpm/src/stdlib_stats_corr.f90. It is listed below.

    module function corr_mask_2_rsp_rsp(x, dim, mask) result(res)
      real(sp), intent(in) :: x(:, :)
      integer, intent(in) :: dim
      logical, intent(in) :: mask(:,:)
      real(sp) :: res(merge(size(x, 1), size(x, 2), mask = 1<dim)&
                          , merge(size(x, 1), size(x, 2), mask = 1<dim))

      integer :: i, j
      real(sp) :: centeri_(merge(size(x, 2), size(x, 1), mask = 1<dim))
      real(sp) :: centerj_(merge(size(x, 2), size(x, 1), mask = 1<dim))
      logical :: mask_(merge(size(x, 2), size(x, 1), mask = 1<dim))

      select case(dim)
        case(1)
          do i = 1, size(res, 2)
            do j = 1, size(res, 1)
             mask_ = merge(.true., .false., mask(:, i) .and. mask(:, j))
             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
                0._sp,&
                mask_)
             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
                0._sp,&
                mask_)

              res(j, i) = dot_product( centerj_, centeri_)&
               /sqrt(dot_product( centeri_, centeri_)*&
                     dot_product( centerj_, centerj_))

            end do
          end do
        case(2)
          do i = 1, size(res, 2)
            do j = 1, size(res, 1)
             mask_ = merge(.true., .false., mask(i, :) .and. mask(j, :))
             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
                0._sp,&
                mask_)
             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
                0._sp,&
                mask_)

              res(j, i) = dot_product( centeri_, centerj_)&
               /sqrt(dot_product( centeri_, centeri_)*&
                     dot_product( centerj_, centerj_))
            end do
          end do
        case default
          call error_stop("ERROR (corr): wrong dimension")
      end select

    end function corr_mask_2_rsp_rsp

The spacing is poor. Why is there a continuation line such as

                0._sp,&
                mask_)

Continuation lines should only be used when a line is too long. Another problem is that merge is used when one of the tsource and fsource arguments require computation. Since Fortran does not mandate short-circuiting, an if block should be used instead. It also looks like the row and column means are computed more often than needed. If they are needed, they should be computed once and stored. The correlation equals the covariance divided by the product of the standard deviations. The standard deviations of the rows or columns should be computed once and stored. Instead, in a line such as

              res(j, i) = dot_product( centeri_, centerj_)&
               /sqrt(dot_product( centeri_, centeri_)*&
                     dot_product( centerj_, centerj_))

the standard deviations for the denominator are computed repeatedly instead of being stored and reused.

Finally,

mask_ = merge(.true., .false., mask(:, i) .and. mask(:, j))

should be written

mask_ =mask(:, i) .and. mask(:, j)

as I wrote in another issue #953.
@Beliavsky
Copy link
Author

Here are other lines in stdlib that do a calculation (computing a mean) in setting tsource or fsource of merge.

The Windows CMD command

c:\fortran\public_domain\github\stdlib>findstr -i -s merge.*mean *.f90 | findstr -v test | findstr -v example

gave the following lines. Presumably could be fixed by changing a few lines of .fypp files.

src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_corr.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:             centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(:, i) - mean(x(:, i), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(:, j) - mean(x(:, j), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centeri_ = merge( x(i, :) - mean(x(i, :), mask = mask_),&
src\temp\stdlib_stats_cov.f90:              centerj_ = merge( x(j, :) - mean(x(j, :), mask = mask_),&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(i, :, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, i, :, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, i, :) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge( (x(:, :, :, i) - mean_)**order,&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(i, :, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, i, :, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, i, :), dp) - mean_)&
src\temp\stdlib_stats_moment_mask.f90:                res = res + merge((real(x(:, :, :, i), dp) - mean_)&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(i, :, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( (x(:, :, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(i, :, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, i, :, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, i, :) - mean)**2,&
src\temp\stdlib_stats_var.f90:                res = res + merge( abs(x(:, :, :, i) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(i, :, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, i, :, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, i, :), dp) - mean)**2,&
src\temp\stdlib_stats_var.f90:              res = res + merge((real(x(:, :, :, i), dp) - mean)**2,&

@jvdp1
Copy link
Member

jvdp1 commented Mar 18, 2025

Thank you @Beliavsky for reviewing the code. Please consider that these parts of stdlib were some of the first added to stdlib, and therefore served as trials and errors for usage of fypp, documentation, ..... I agree that some parts of these codes would not pass the current review standard of stdlib. Please also consider you are reviewing the preprocessed files, and not the fypp files. This might explained the observed repetitions of code, or strange formatting.

The spacing is poor. Why is there a continuation line such as

            0._sp,&
            mask_)

Continuation lines should only be used when a line is too long.

This is the result of the preprocessing. the formatting makes sense in the fypp files (at least to me).

Another problem is that merge is used when one of the tsource and fsource arguments require computation. Since Fortran does not mandate short-circuiting, an if block should be used instead.

if you refer to merge( x(:, j) - mean(x(:, j), mask = mask_),& 0._sp,& mask_)
my experience is that some compilers can better optimize this approach than a if block. It is also more consice and easier to read than if blocks and loops. where intrinsics could have been used.

It also looks like the row and column means are computed more often than needed. If they are needed, they should be computed once and stored. The correlation equals the covariance divided by the product of the standard deviations. The standard deviations of the rows or columns should be computed once and stored. Instead, in a line such as

          res(j, i) = dot_product( centeri_, centerj_)&
           /sqrt(dot_product( centeri_, centeri_)*&
                 dot_product( centerj_, centerj_))

the standard deviations for the denominator are computed repeatedly instead of being stored and reused.

In this case, I don't think it could be done as you proposed because the sd will depend on each combinasion mask(:, i) .and. mask(:, j))

mask_ = merge(.true., .false., mask(:, i) .and. mask(:, j))

should be written

mask_ =mask(:, i) .and. mask(:, j)

You are right, it should be fixed. Could you open a PR, please?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
2 participants
@jvdp1 @Beliavsky