裏 RjpWiki

Julia ときどき R, Python によるコンピュータプログラム,コンピュータ・サイエンス,統計学

Julia のベクトル(その 2)

2021年02月08日 | ブログラミング

並べ替え

d = [3, 2, 5, 1, 4];
sort(d) # 1, 2, 3, 4, 5; 昇順ソート
sort(d, rev=true) # 5, 4, 3, 2, 1; 降順ソート

値のある場所のインデックス

一番大きい要素のインデックス,二番目に大きい要素のインデックス,...

sortperm(d) # 4, 2, 5, 3, 1; 昇順ソートするときに,順に取り出す要素の位置
d[sortperm(d)] # sortperm(d) を使ってソートするのはこういうこと
sortperm(d, rev=true) # 3, 5, 1, 2, 4; 降順ソートするときに,順に取り出す要素の位置
d[sortperm(d, rev=true)] # 5, 4, 3, 2, 1; sortperm(d) を使ってソートするのはこういうこと
sort!(d) # 1, 2, 3, 4, 5; 破壊的ソート。d がインプレースでソートされるので,代入しなくてよい。
d # 1, 2, 3, 4, 5; 内容が変わっている

統計関数

using Statistics
using Random;
Random.seed!(12345); # 乱数の種
x = randn(100);

最大値

maximum(x) # 2.1307390315181536; 最大値
argmax(x)  # 55; 最大値のあるインデックス
findmax(x) # (2.1307390315181536, 55)); 最大値とそのインデックス

最小値

minimum(x) # -2.3280152457586003; 最小値
argmin(x)  # 73; 最小値のあるインデックス
findmin(x) # (-2.3280152457586003, 73); 最小値とそのインデックス

最大値と最小値

maximu(), minimu() 両方やるよりも遅い

extrema(x) # (-2.3280152457586003, 2.1307390315181536);

sum(x) # 0.9821131440360702

平均値

mean(x) # 0.009821131440360701

不偏分散,分散

var(x) # 0.9954011110241598
var(x, corrected=false) # 0.9854470999139182

不偏標準編纂,標準偏差

std(x) # 0.9976979056929807
std(x, corrected=false) # 0.9926968821921011

中央値

median(x) # 0.13578110308214625
median!(x) # 0.13578110308214625; 破壊的(インプレース)関数

要約統計量

describe(x)
# Summary Stats:
# Length:         100
# Missing Count:  0
# Mean:           0.009821
# Minimum:        -2.328015
# 1st Quartile:   -0.684464
# Median:         0.135781
# 3rd Quartile:   0.690881
# Maximum:        2.130739
# Type:           Float64

y = [2, 4, 1, 3, 5];

diff(y) # 2, -3, 2, 2

累和

cumsum(y) # 2, 6, 7, 10, 15
[sum(y[1:i]) for i = 1:length(y)] # 2, 6, 7, 10, 15
z = similar(y);
cumsum!(z, y); # 2, 6, 7, 10, 15; 破壊的(インプレース)関数
y; # 2, 4, 1, 3, 5
z; 2, 6, 7, 10, 15

y = [2, 4, 1, 3, 5];
prod(y) # 120

cumprod, cummin, cummax

[prod(y[1:i]) for i = 1:length(y)] # 2, 8, 8, 24, 120; cumprod() に相当
[minimum(y[1:i]) for i = 1:length(y)] # 2, 2, 1, 1, 1; cummin() に相当
[maximum(y[1:i]) for i = 1:length(y)] # 2, 4, 4, 4, 5; cummax() に相当

不偏共分散,共分散

a = [1, 3, 4, 2, 6];
b = [3, 2, 5, 4, 3];
cov(a, b) # 0.15000000000000002; 不偏共分散
cov(a, b, corrected=false) # 0.12000000000000002; 共分散

ピアソンの積率相関係数

cor(a, b) # 0.06839411288813299; ピアソンの積率相関係数
A = randn(100, 3);
cor(A) # ピアソンの積率相関係数行列
# 3×3 Array{Float64,2}:
#  1.0         0.0227498   0.18433
#  0.0227498   1.0        -0.0643971
#  0.18433    -0.0643971   1.0

ノルム

using LinearAlgebra
norm(a) # 8.12403840463596; ノルム

内積

dot(a, b) # 55; 内積
sum(a .* b) # 55; 内積

外積

a = [1, 2, 3];
b = [10, 20, 30, 40];
a * b' # 外積
# 3×4 Array{Int64,2}:
#  10  20  30   40
#  20  40  60   80
#  30  60  90  120

using StatsBase  # geomean, harmmean, quantile, corspearman, corkendall

分位数(パーセンタイル値)

quantile(x)
  # -2.3280152457586003; 最小値
  #  -0.6844644232186423; 25%タイル値
  #   0.13578110308214625; 中央値(50%タイル値)
  #   0.690881410484019; 75%タイル値
  #   2.1307390315181536; 最大値
quantile(x, 0.00) # -2.3280152457586003; 最小値
quantile(x, 0.25) # -0.6844644232186423; 25%タイル値
quantile(x, 0.50) # 0.13578110308214625; 中央値(50%タイル値)
quantile(x, 0.75) # 0.690881410484019; 75%タイル値
quantile(x, 1.00) # 2.1307390315181536; 最大値
quantile!(x, 0.5) # 0.13578110308214625; 破壊的(インプレース)関数

幾何平均

geomean([2, 1, 3, 2, 1]) # 1.6437518295172258; 幾何平均

調和平均

harmmean([2, 1, 3, 2, 1]) # 1.5; 調和平均

スピアマンの順位相関係数

a = [1, 3, 4, 2, 6];
b = [3, 2, 5, 4, 3];
corspearman(a, b) # 0.10259783520851543; スピアマンの順位相関係数
A = randn(10, 3);
corspearman(A); # スピアマンの順位相関係数行列
# 3×3 Array{Float64,2}:
#   1.0         -0.430303   0.00606061
#  -0.430303     1.0       -0.139394
#   0.00606061  -0.139394   1.0

ケンドールの順位相関係数

corkendall(a, b) # 0.10540925533894598; ケンドールの順位相関係数
corkendall(A); # ケンドールの順位相関係数行列
# 3×3 Array{Float64,2}:
#   1.0        -0.288889    0.0222222
#  -0.288889    1.0        -0.0222222
#   0.0222222  -0.0222222   1.0

集合

v = [2, 2, 4, 4, 4, 4, 1, 3, 3, 3, 5, 5, 5, 5, 5];
A = Set(v); # 4, 2, 3, 5, 1; 集合(要素の順序はない。重複もない。)
B = Set([1, 3, 5, 7, 9]); 7, 9, 3, 5, 1
union(A, B) # 7, 4, 9, 2, 3, 5, 1; 和集合
intersect(A, B) # 3, 5, 1; 積集合
setdiff(A, B) # 4, 2; 差集合
symdiff(A, B) # 7, 4, 9, 2; 対称差集合
A == B # false; 等しい
1 in A # true; 集合に含まれる
0 in A # false
Set() # (); 空集合
empty(A) # (); 空集合にする

要素の置換

v = [2, 4, 1, 3, 5]; # 2, 4, 1, 3, 5
v[[1, 3, 5]] = [10, 30, 50];
v; # 10, 4, 30, 3, 50

ベクトルの結合

u = [1, 2, 3];
v = [4, 5];
w = [6, 7, 8, 9];
vcat(u, v, w); # 1, 2, 3, 4, 5, 6, 7, 8, 9
append!(u, v); # インプレースで結合
append!(u, w);
u; # 1, 2, 3, 4, 5, 6, 7, 8, 9

要素・ベクトルの挿入

a = [1, 2, 3];
pushfirst!(a, 0) # 0, 1, 2, 3; インプレースで先頭に挿入

a = [1, 2, 3];
push!(a, 99) # 1, 2, 3, 99; インプレースで末尾に挿入(追加)

a = [1, 2, 3];
insert!(a, 2, 99) # 1, 99, 2, 3; インプレースでインデックス位置に挿入

u = [1, 2, 3];
v = [4, 5];
[insert!(u, 3, v2) for v2 in reverse(v)]
u; # 1, 2, 4, 5, 3

ベクトルの要素の削除

v = [1, 2, 3, 4, 5];
popfirst!(v) # 1; インプレースで先頭の要素を削除して,返す
v; # 2, 3, 4, 5

v = [1, 2, 3, 4, 5];
pop!(v) # 5; インプレースで最後の要素を削除して,返す
v; # 1, 2, 3, 4

v = [1, 2, 3, 4, 5];
deleteat!(v, 3); # インプレースでインデックス位置の要素を削除(delete at)
v; # 1, 2, 4, 5

empty!(v); # 空ベクトルにする
v; # Int64[]

要素ごとの演算

v = [2, 4, 1, 3, 5];
w = [1, 3, 5, 7, 9];
v + w;  # 3, 7, 6, 10, 14
v .+ w; # 3, 7, 6, 10, 14
v - w;  # 1, 1, -4, -4, -4
v .- w; # 1, 1, -4, -4, -4
v .* w; # 2, 12, 5, 21, 45
v * w; # エラー MethodError: no method matching *(::Array{Int64,1}, ::Array{Int64,1})
v ./ w; # 2.0, 1.3333333333333333, 0.2, 0.42857142857142855, 0.5555555555555556
v / w # 注! 5 x 5 配列
# 5×5 Array{Float64,2}:
#  0.0121212   0.0363636  0.0606061  0.0848485  0.109091
#  0.0242424   0.0727273  0.121212   0.169697   0.218182
#  0.00606061  0.0181818  0.030303   0.0424242  0.0545455
#  0.0181818   0.0545455  0.0909091  0.127273   0.163636
#  0.030303    0.0909091  0.151515   0.212121   0.272727
v .+ 1; # 3, 5, 2, 4, 6
v + 1; # エラー MethodError: no method matching +(::Array{Int64,1}, ::Int64)
v * 2; # 4, 8, 2, 6, 10
v .* 2; # 4, 8, 2, 6, 10
v .^ 2; # 4, 16, 1, 9, 25
v ^ 2; # エラー MethodError: no method matching ^(::Array{Int64,1}, ::Int64)
1 ./ v # 0.5, 0.25, 1.0, 0.333, 0.2
1 / v # 0.0363636  0.0727273  0.0181818  0.0545455  0.0909091 ?? どういう計算?
v .< 3 # 1, 0, 1, 0, 0
mod.(v, 2) # 0, 0, 1, 1, 1
mod.(v, 2) .== 0; # 1, 1, 0, 0, 0
sqrt.(v); # 1.41421, 2.0, 1.0, 1.73205, 2.23607
exp.(v); # 7.38906, 54.5982, 2.71828, 20.0855, 148.413

定数

円周率

pi # π = 3.1415926535897... 円周率
π # π = 3.1415926535897...; # \pi<tab> で入力
typeof(π) # Irrational{:π}
1pi # 3.141592653589793
typeof(1pi) # Float64
BigFloat(pi) # 3.141592653589793238462643383279502884197169399375105820974944592307816406286198
BigFloat(1) + pi # 4.141592653589793238462643383279502884197169399375105820974944592307816406286198
BigFloat(1) + Float64(pi) # 4.141592653589793115997963468544185161590576171875
1 + pi # 4.141592653589793
Float64(pi) # 3.141592653589793

1024 ビットにする

setprecision(BigFloat, 1024) do
    BigFloat(pi)
end
# 3.14159265358979323846264338...24586997

ネイピア数

# ℯ = 2.7182818284590...; \euler<tab> で入力; ネイピア数
typeof(ℯ) # Irrational{:ℯ}
BigFloat(ℯ) # 2.718281828459045235360287471352662497757247093699959574966967627724076630353555
1ℯ # 2.718281828459045

 

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

Julia のベクトル(その 1)

2021年02月08日 | ブログラミング

ベクトル

Julia では,配列(特に1次元のときはベクトル,2次元のときは行列とも呼ぶ)は [ ] で囲んで使用する。

v = [1, 2, 3, 4, 5];
v # 1, 2, 3, 4, 5

対話モードで実行する場合(RERL; Run Evaluate Print Loop と呼ばれる),行の最後に ';' が付いているとその行の評価結果を表示しない。C 言語などのように,行末に必ず必要な ';' ではない。

typeof(v) # Array{Int64,1}

typeof
関数は,オブジェクトの型を表示する。上のように小数点なしの数を ',' で区切って定義すると整数ベクトルになる。Array は配列,Int64 は要素が整数,1 は次元数が 1 すなわちベクトルであることを意味する。

ちなみに,配列の次元サイズは size() で求めることができる。

size(v) # (5,)

次元が 1 のときは,上の例のように示される。
最後の ',' は (5) だと数値の 5 がカッコでくくられたもの(数式)と区別がつかないのでこのように示される。
つまり,この場合の '( )' は数式の中での '( )' ではなく,これ全体がタプル(複数の要素のまとまり)という型であることを意味している。

要素に一つでも小数点がついていると実数ベクトルになる。

v2 = [1.0, 2, 3, 4, 5]
; # 1.0, 2.0, 3.0, 4.0, 5.0

整数ベクトルを実数ベクトルにするときは,以下のように型を指定して変換することができる。

v3 = Array{Float64,1}(v); # 1.0, 2.0, 3.0, 4.0, 5.0

Float64 で,実数型であることを示している。

typeof(v3) # Array{Float64,1}

または,

v4 = convert.(Float64, v); # 1.0, 2.0, 3.0, 4.0, 5.0

convert は最初の引数で示した型に変換する関数。'(' の前の '.' は,第 2 引数が複数の要素を持つ場合,それぞれの要素に作用することを指示する。もし,この '.' がないとエラーになる。

空白で区切ると横ベクトルになる。

v5 = [10 20 30 40 50];
# 1×5 Array{Int64,2}:
#  10  20  30  40  50

しかしこれは 1 行 n 列の配列(行列)である。typeof() で型を見ると次元数が 2 になっていることで確認できる。

typeof(v5) # Array{Int64,2}
size(v5)   # (1, 5)

v の定義のときのように ',' で区切って定義した縦ベクトルを n 行 1 列の配列(行列)にするには,配列の形を変える関数 reshape を使って次のようにする。

reshape(v, 5, 1) # 1, 2, 3, 4, 5

第 2,第 3 引数の 5, 1 は 結果として返す配列の第 1,第 2 次元の大きさを示すので,5 行 1 列の配列になる。
一般的には以下のようにするとよい。length(v)v の長さ,':' は 2 次元配列の場合は,1 次元目の大きさが決まれば自動的に決まるので,特に示す必要がないときに指定するとよい。

reshape(v, length(v), :)

逆に ' ' で区切った 1 行 n 列の配列をベクトルにしたいときには vec 関数を使う。

vec([10 20 30 40 50]) # 10, 20, 30, 40, 50

要素が文字列である配列(ベクトル,行列)も同じように定義できる。typeof 関数で示される String は要素が文字型であることを示している。

v8 = ["A", "B", "C", "D", "E"]; # "A", "B", "C", "D", "E"
typeof(v8) # Array{String,1}

Julia では,文字列は長さが 1 であっても,二重引用符で囲まなければならない。
以下のように,1 文字を一重引用符で囲むと文字列型ではなく文字型 Char になる。逆に言えば,長さが2以上の文字列を一重引用符で囲むとエラーになる。

v7 = ['1', 'A', 'b'];
# 3-element Array{Char,1}:
#  '1': ASCII/Unicode U+0031 (category Nd: Number, decimal digit)
#  'A': ASCII/Unicode U+0041 (category Lu: Letter, uppercase)
#  'b': ASCII/Unicode U+0062 (category Ll: Letter, lowercase)
typeof(v7) # Array{Char,1}

規則性のあるベクトル

等差数列

collect(1:1:10) # 1 〜 10,初項 1,公差 1; 第2引数が公差
collect(1:10) # 1 〜 10,初項 1,公差 1; 第2引数が1の場合は省略可
collect(10:-1:1) # 10 〜 1,初項 10,公差 -1; 降順の場合も公差は明示する
collect(1:5:20) # 1, 6, 11, 16,初項 1,公差 5; 第2引数が公差
collect(0:0.2:3) # 初項 0, 公差 0.2の数列; 最終項は 3

collect は要素を実際にメモリー上に確保する。
for i = 1:10 のような使い方(イテレータ)のときには,メモリー上に確保しない。

range(1, 10, step=2) # 1:2:9; イテレータを作成する
range(1, 10, length=20) # 1.0:0.47368421052631576:10.0
collect(range(1, 10, length=3)); # 1.0, 5.5, 10.0

繰り返し

repeat([0], 5) # 0 が 5個; repeat(0, 5) はエラーになる
repeat([0, 1, 2], 5) # 0, 1, 2 が 5 回繰り返されるベクトル
repeat(["foo", "bar", "baz"], 5) # "foo", "bar", "baz" が 5 回繰り返されるベクトル

A = repeat([10 20 30], inner=(1,2,1)) # 10  10  20  20  30  30
B = repeat([10 20 30], inner=(1,4,1)) # 10  10  10  10  20  20  20  20  30  30  30  30

inner を指定する repeat 関数,repeat(v, inner=(1, m, 1))1 × (length(v)×m) × 1 の3次元配列である。

size(A) # (1, 6, 1)
size(B) # (1, 12, 1)

vec 関数で,ベクトルに変換できる。

vec(A) # 10  10  20  20  30  30
vec(B) # 10  10  10  10  20  20  20  20  30  30  30  30

ゼロベクトル

zeros(5) # 実数 0.0 が 5 個のベクトル
zeros(Int64, 5) # 整数 0 が 5 個; 64 ビット環境では zeros(Int, 5) と同じ

1 ベクトル

ones(5) # 実数 1.0 が 5 個のベクトル; fill(1.0, 5) と同じ
ones(Int64, 5) # 整数 1 が 5 個; 64 ビット環境では zeros(Int, 5) と同じ

任意の要素の繰り返しを持つベクトル

fill(0, 5) # zeros(Int64, 5) と同じ
fill(1, 5) # ones(Int64, 5) と同じ
fill(0.0, 5) # zeros(5) と同じ
fill(1.0, 5) # ones(1) と同じ
fill("a", 5) # "a" が 5 個のベクトル

乱数ベクトル

using Random;
Random.seed!(12345); # 乱数の種の設定。本当に乱数が必要なときは使わない。

以下では,途中で乱数を使わないときには,示されたのと同じ結果になる。

無作為抽出

ベクトルから要素をランダムに取り出す。

rand(1:5, 10) # 3, 4, 1, 2, 2, 5, 5, 4, 5, 2; 要素 1〜5 から 10 個を無作為抽出したベクトル
rand(["foo", "bar", "baz"], 5) # "baz", "bar", "bar", "foo", "bar"

乱数

一様乱数を要素とする,長さ 3 のベクトル

rand(3) # 0.36007770681119133, 0.41511087872028374, 0.8325976767993883

正規乱数を要素とする,長さ  3 のベクトル

randn(3) # -0.55971756622102, -1.4990404969589635, 1.6739451586233716

要素の取得

Julia では,R と同様先頭要素のインデックスは 1 である。

v11 = [10, 20, 30, 40, 100]; # 10, 20, 30, 40, 100
v11[1] # 10
v11[begin] # 10; 'begin' はインデックス 1 を表す
v11[end] # 100; 'end' は最後のインデックスを表す
v11[end-2] # 30; 'begin', 'end' を算術演算してインデックスにできる
v11[[2,3,end]] # 20, 30, 100 # ベクトルで指定して要素を取り出すことができる
v11[[1,2,1,3,5,5]] # 10, 20, 10, 30, 100, 100; 何回でもどこからでも取り出せる
v11[1:2:end] # 10, 30, 100; イテレータで要素を取り出すことができる
v11[v11 .> 15] # 20, 30, 40, 100; 条件式(比較演算)を満たす要素を取り出せる

'.>''.' は複数の要素を持つオブジェクトの各要素に適用するために必要である。
なければエラーになる。

BitArray 型のベクトル

v11 .> 15 # BitArray 型のベクトル [0, 1, 1, 1, 1]

配列のインデックスに指定すると,1 に対応するオブジェクトの要素が取り出される。

v11[15 .< v11 .< 35] # 20, 30; 条件式(論理演算))を満たす要素を取り出せる

'15 .< v11 .< 35' '15 .< v11 .& v11 .< 35' と同じではない。

以下の 2 つは 同じことを表していそうに思うが,結果が違う。

v11[15 .< v11 .& v11 .< 35] # 20, 30; 両方の条件を満たす要素を取り出す。 & は and
v11[v11 .> 15 .& v11 .< 35] # 20, 30, 40, 100; 両方の条件を満たす要素を取り出す。 & は and

'15 .< v11 .< 35''15 .< v11 .& v11 .< 35' と同じになる。'.&' とすることに注意。
'15 .< v11''v11 .> 15' は同じ意味なので,
'v11[15 .< v11 .& v11 .< 35]''v11[v11 .> 15 .& v11 .< 35]' の結果が違うのは一見納得ができない。

実は,ここでの '&' はビット演算子であり,比較演算子 '.<''.>' より優先順位が高い。
'v11 .& v11''[10, 20 ,30, 40, 100]',で
'v11[15 .< [10, 20 ,30, 40, 100] .< 35]''[20, 30]' になる。

一方,'15 .& v11''[10, 4, 14, 8, 4]' で,
'v11[v11 .> [10, 4, 14, 8, 4] .< 35]''20, 30, 40, 100' になる。当然,予期せぬ間違った答えである。

以下のように,比較演算を先に行うようにカッコを補えば,正しい答えになる。

v11[(v11 .> 15) .& (v11 .< 35)] # 20, 30

また,'v11[15 .< v11 .& v11 .< 35]' もプログラムを書いた意図からは 'v11[(15 .< v11) .& (v11 .< 35)]' とすべきである。

以下のような単純な論理式の例だとわかりやすいかもしれない。

1 < 3 & 4 < 5 # false
(1 < 3) & (4 < 5) # true

1 > 3 | 5 < 6 # false
(1 > 3) | (5 < 6) # true

Python の場合も R の場合も,論理演算子のが出現する可能性のあるところにビット演算子も出現する可能性はないので,Julia のときのような誤解は生じないのである。

R にはビット演算はないので,'&' は論理演算子として解釈されるので,比較演算をカッコでくくらなくても正しい答えになる。

using RCall
R"""
v11 = c(10, 20, 30, 40, 100)
v11[v11 >= 15 & v11 <= 35] # 20, 30
"""

要素への代入

インデックスを指定してスカラー,ベクトルを代入することができる。

v = collect(1.:10);
v[1] = 10;
v[2:5] = [20, 30, 40 ,50];
v[9] = NaN;
v # 10.0  20.0  30.0  40.0  50.0  6.0  7.0  8.0  NaN  10.0

脇道

欠損値は 'missing' で表すが,上のようなベクトル v の要素として代入はできない。

'v[8] = missing' はエラーになる。

'missing' も許す型(たとえば 'Union{Missing, Int64}')のオブジェクトならば,可能である。

x = Array{Union{Missing, Int64},1}([3,2,1,4,3,5,missing,6]);
x[4] = missing;
x #  3, 2, 1, missing, 3, 5, missing, 6

たとえば,以下のようにすると 'missing' を初期値として持つ長さ 10 のベクトルを定義できる。

y = Array{Union{Missing, Float64},1}(undef, 10)

1, 4, 7, 10 番目の要素を二乗する。

for i =1:3:length(y)
  y[i] = i^2
end
y; # 1, 2, 1, 16, 5

'missing' を取り除いたベクトルを求める。

z = collect(skipmissing(y)); # 1, 4, 1, 16, 5
using Statistics # 単に mean を使うだけなのにこれが必要
mean(z) # 41.5; missing のないベクトルの平均値
mean(y) # missing
mean(skipmissing(y)) # 41.5; missing を除いたベクトルの平均値

Julia に NA はない。

NA; # UndefVarError: NA not defined

Julia の NaN は?

NaN # NaN
NaN == missing # missing
NaN * 10.3 # NaN
NaN * Inf # NaN
NaN & missing # MethodError: no method matching &(::Float64, ::Missing)

ベクトルの属性

a = [1, 2, 3, 4, 5];

長さ

length(a) # 5

検索

要素はあるか?

2 in a # true
9 in a # false

何個あるか?

b = [0, 1, 2, 1, 2, 3, 4, 3, 5, 2, 2, 1, 1, 5, 0];

count の第一引数は関数である。
比較演算式,論理演算式も関数なので,count の第 1 引数になれる。

5 < 10 # true
<(5, 10) # true

count(isequal(1), b) # 4; 1 と等しいか?
count(b .== 1) # 4; b は 1 と等しいか?
count(isequal.(b, 1)) # 4; b と 1 が等しいか? 等しければ 1 等しくなければ 0
count(iseven, b) # 7; 偶数か?
count(mod.(b, 2) .== 0) # 7; 偶数か?
count(isless.(b, 2)) # 6; b は 2 より小さいか? isgreater() はない
count(b .< 2) # 6; b は 2 より小さいか?

最大値,最小値は何か?どこにあるか?

argmax(b) # 9; 最大値が最初に現れるインデックス
b[argmax(b)] # 5; 最大値
maximum(b) # 5; 最大値
argmin(b) # 1; 最小値が最初に現れるインデックス
b[argmin(b)] # 0; 最小値
minimum(b) # 0; 最小値

ある値がどこにあるか?

findall 関数の第 1 引数は,関数である。論理演算子も関数である。λ関数(無名関数)でもよい。

findall(==(1), b); # 2, 4, 12, 13; 1 と等しい要素のインデックス
findall(==(b[argmax(b)]), b) # 9, 14; 最大値が現れる全てのインデックス
findall(==(minimum(b)), b) # 9, 14; 最大値が現れる全てのインデックス
findall(==(b[argmin(b)]), b) # 1, 15; 最小値が現れる全てのインデックス
findall(==(minimum(b)), b) # 1, 15; 最小値が現れる全てのインデックス

findall(isodd, b) # 2, 4, 6, 8, 9, 12, 13, 14; 要素が奇数のインデックス
b[findall(isodd, b)] # 1, 1, 3, 3, 5, 1, 1, 5; その要素の列挙
findall(x -> 2 < x < 4, b) # 6, 8; 要素を x として 2 < x < 4 を満たす要素のインデックス(λ関数)
b[findall(x -> 2 < x < 4, b)] # 3, 3
c = findall(x -> sqrt(x) == floor(Int, sqrt(x)), b) # 1, 2, 4, 7, 12, 13, 15; 要素が平方数であるインデックス
b[c] # 0, 1, 1, 4, 1, 1, 0; その要素

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

Julia で注意しておかないといけないこと(その 1)

2021年02月08日 | ブログラミング

以下の 2 つの論理式は true にはならない

1 < 3 & 4 < 5 # false
1 > 3 | 5 < 6 # false

理由を知りたい人は,ずっと下にスクロールする

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

以下のように書かなければならない

(1 < 3) & (4 < 5) # true
(1 > 3) | (5 < 6) # true

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

整数平方根 isqrt

2021年02月08日 | ブログラミング

整数 n に対して a^2 <= n を満たす最大の整数 a を求めるときは,
int(sqrt(n)) ではなくisqrt を使う

ということであるが,math に isqrt があるのは知らなくて,自分で関数を書いたりした。ニュートン法なのでそこそこの速度はあったのだけど。

from math import isqrt
from time import time

def f(n):
    a = isqrt(n)
    return a**2 <= n < (a+1)**2

start = time(); print(f(2**12345678), time()-start)
# True 37.5530788898468

Julia にも isqrt はあるが,実引数 n の指定には BigInt が使われるように定義すること。まあ,それだけだから,簡単といえば簡単だ。

function f(n)
    a = isqrt(n)
    a^2 <= n < (a+1)^2
end

@time f(BigInt(2)^12345678)
# 0.178914 seconds (1.44 k allocations: 51.893 MiB, 0.52% gc time)
# true

Python より 200 倍速い。

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

PVアクセスランキング にほんブログ村

PVアクセスランキング にほんブログ村