小编JEq*_*hua的帖子

所以我发现了这个:

转换MATLAB代码时,可能需要首先将矩阵重新整形为线性序列,执行一些索引操作,然后重新整形.由于重塑(通常)会在同一存储上生成视图,因此应该可以相当有效地执行此操作.

请注意,Numpy中reshape使用的扫描顺序默认为'C'顺序,而MATLAB使用Fortran顺序.如果您只是简单地转换为线性序列,那么这无关紧要.但是如果要从依赖于扫描顺序的MATLAB代码转换重构,那么这个MATLAB代码:

z = reshape(x,3,4);

Run Code Online (Sandbox Code Playgroud)

应该成为

z = x.reshape(3,4,order='F').copy()

Run Code Online (Sandbox Code Playgroud)

在Numpy.

mafs当我在MATLAB中执行时,我有一个多维16*2数组:

mafs2 = reshape(mafs,[4,4,2]) 

我得到的东西不同于我在python中做的事情:

mafs2 = reshape(mafs,(4,4,2))

甚至

mafs2 = mafs.reshape((4,4,2),order='F').copy()

对此有何帮助？谢谢你们.

我无法准确理解反射模式如何处理我的数组.我有这个非常简单的数组:

import numpy as np
from scipy.ndimage.filters import uniform_filter
from scipy.ndimage.filters import median_filter

vector = np.array([[1.0,1.0,1.0,1.0,1.0],[2.0,2.0,2.0,2.0,2.0],[4.0,4.0,4.0,4.0,4.0],[5.0,5.0,5.0,5.0,5.0]])

print(vector)

[[1. 1. 1. 1. 1.] [2. 2. 2. 2. 2.] [4. 4. 4. 4. 4.] [5. 5. 5. 5. 5.]]

应用窗口大小为3的均匀(平均)过滤器,我得到以下结果:

filtered = uniform_filter(vector, 3, mode='reflect')

print(filtered)

[[1.33333333 1.33333333 1.33333333 1.33333333 1.33333333] [2.33333333 2.33333333 2.33333333 2.33333333 2.33333333] [3.66666667 3.66666667 3.66666667 3.66666667 3.66666667] [4.66666667 4.66666667 4.66666667 4.66666667 4.66666667]]

如果我尝试手动复制练习,我可以得到这个结果.原始矩阵为绿色,窗口为橙色,结果为黄色.白色是"反映"的观察结果.

结果是:

但是当我尝试4或5的窗口大小时,我无法复制结果.

filtered = uniform_filter(vector, 4, mode='reflect')

print(filtered)

[[1.5 1.5 1.5 1.5 1.5] [2. 2. …

我一直在玩numba和numexpr试图加速一个简单的元素矩阵乘法.我无法获得更好的结果,它们基本上(速度方向)等同于numpys乘法函数.这个地区有人有运气吗？我使用numba和numexpr是错误的(我对此很新)或者这是一个不好的方法来尝试加快速度.这是一个可重现的代码,谢谢你的高级:

import numpy as np
from numba import autojit
import numexpr as ne

a=np.random.rand(10,5000000)

# numpy
multiplication1 = np.multiply(a,a)

# numba
def multiplix(X,Y):
    M = X.shape[0]
    N = X.shape[1]
    D = np.empty((M, N), dtype=np.float)
    for i in range(M):
        for j in range(N):
            D[i,j] = X[i, j] * Y[i, j]
    return D

mul = autojit(multiplix)
multiplication2 = mul(a,a)

# numexpr
def numexprmult(X,Y):
    M = X.shape[0]
    N = X.shape[1]
    return ne.evaluate("X * Y")

multiplication3 = numexprmult(a,a) 

我有一个5000*5000numpy数组,我想要计算大小为25的窗口的Kurtosis.我尝试将scipys自己的kurtosis功能放入generic_filter发现中,ndimage.filters如下所示:

import numpy as np

from scipy.stats import kurtosis
from scipy.ndimage.filters import generic_filter

mat = np.random.random_sample((5000, 5000))

kurtosis_filter = generic_filter(mat, kurtosis, size=25, mode='reflect') 

这永远不会结束,我不确定它会给出正确的答案.所以我的第一个问题是,这是否是使用generic_filterscipy函数的正确方法.如果它恰好是正确的,那么它对我来说太慢了.所以我的下一个问题是,如果有更快的方法来实现这一目标？例如,考虑标准偏差,您可以简单地执行以下操作:

usual_mean = uniform_filter(mat, size=25, mode='reflect')
mean_of_squared = uniform_filter(np.multiply(mat,mat), size=25, mode='reflect')
standard_deviation = (mean_of_squared - np.multiply(usual_mean,usual_mean))**.5

这是非常快速的,简单来自$\sigma ^ 2 = E [(X - \mu)^ 2] = E [X ^ 2] - (E [X])^ 2 $的事实.

我一直在使用MASS包,可以使用image和par(new = TRUE)绘制两个双变量法线,例如:

# lets first simulate a bivariate normal sample
library(MASS)
bivn <- mvrnorm(1000, mu = c(0, 0), Sigma = matrix(c(1, .5, .5, 1), 2))
bivn2 <- mvrnorm(1000, mu = c(0, 0), Sigma = matrix(c(1.5, 1.5, 1.5, 1.5), 2))

# now we do a kernel density estimate
bivn.kde <- kde2d(bivn[,1], bivn[,2], n = 50)
bivn.kde2 <- kde2d(bivn2[,1], bivn[,2], n = 50)

# fancy perspective
persp(bivn.kde, phi = 45, theta = 30, shade = .1, border = NA)
par(new=TRUE)
persp(bivn.kde2, …

我正在尝试向量化滑动窗口操作.对于1-d案例,一个有用的例子可以遵循:

x= vstack((np.array([range(10)]),np.array([range(10)])))

x[1,:]=np.where((x[0,:]<5)&(x[0,:]>0),x[1,x[0,:]+1],x[1,:])

index <5的每个当前值的n + 1值.但我得到这个错误:

x[1,:]=np.where((x[0,:]<2)&(x[0,:]>0),x[1,x[0,:]+1],x[1,:])
IndexError: index (10) out of range (0<=index<9) in dimension 1

奇怪的是,我不会得到n-1值的这个错误,这意味着索引小于0.它似乎并不介意:

x[1,:]=np.where((x[0,:]<5)&(x[0,:]>0),x[1,x[0,:]-1],x[1,:])

print(x)

[[0 1 2 3 4 5 6 7 8 9]
 [0 0 1 2 3 5 6 7 8 9]]

有没有办法解决？我的方法完全错了吗？任何意见将不胜感激.

编辑:

这就是我想要实现的目标,我将矩阵展平为一个numpy数组,我想要计算每个单元格的6x6邻域的平均值:

matriz = np.array([[1,2,3,4,5],
   [6,5,4,3,2],
   [1,1,2,2,3],
   [3,3,2,2,1],
   [3,2,1,3,2],
   [1,2,3,1,2]])

# matrix to vector
vector2 = ndarray.flatten(matriz)

ncols = int(shape(matriz)[1])
nrows = int(shape(matriz)[0])

vector = np.zeros(nrows*ncols,dtype='float64')


# Interior pixels
if ( (i % ncols) != 0 and …

我有以下矩阵sigma和sigmad:

西格玛:

    1.9958   0.7250
    0.7250   1.3167

sigmad:

    4.8889   1.1944
    1.1944   4.2361

如果我试图解决python中的广义特征值问题,我得到:

    d,V = sc.linalg.eig(matrix(sigmad),matrix(sigma))

五:

    -1     -0.5614
    -0.4352    1

如果我尝试在matlab中解决ge问题,我会得到:

    [V,d]=eig(sigmad,sigma)

五:

    -0.5897    -0.5278
    -0.2564    0.9400

但是d确实很重要.

我试图从R调用winBUGS来估计逻辑回归.我正在使用以下代码:

# Directorio de trabajo
setwd("~/3 Diplomado/7 Bayesiana/8t1")

# paquete para hablarse con WinBUGS desde R
library(R2WinBUGS)

# cargamos datos
reg <- read.table("enf.csv", header = TRUE, sep = ",")
edad <- reg$edad
enfer <- reg$efer
n <- length(reg$edad)

# Primeras filas de los datos
head(reg)

# Nombres de los datos para alimentar al modelo en WinBUGS
datos <- list("edad","enfer", "n")

# Construimos el modelo

modelo <- function(){
  for (i in 1:n) {
    enfer[i] ~ dbin(theta[i], 1)
    logit(theta[i]) < - beta0 …