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

making kmeans take AbstractMatrix instead of Matrix #138

Merged
merged 13 commits into from
Jan 29, 2019
Merged

making kmeans take AbstractMatrix instead of Matrix #138

Changes from 1 commit
Commits
Show all changes
13 commits
Select commit Hold shift + click to select a range
57c172a
making kmeans take AbstractMatrix instead of Matrix, see #136 and #129
tlienart Jan 25, 2019
4c7f218
PR update based on feedback + cosmetics
tlienart Jan 26, 2019
2e95afb
Update src/kmeans.jl
alyst Jan 28, 2019
be3575d
Update src/kmeans.jl
alyst Jan 28, 2019
d3f167a
Update src/kmeans.jl
alyst Jan 28, 2019
a8ea8c7
Update src/kmeans.jl
alyst Jan 28, 2019
c5bda5b
Update src/kmeans.jl
alyst Jan 28, 2019
3f8b15b
Update src/kmeans.jl
alyst Jan 28, 2019
0d272e9
applying the changes following PR reviews
tlienart Jan 28, 2019
9cb02f7
reverting a few changes that seemed to cause spurious allocations
tlienart Jan 28, 2019
16ae803
type changes, tying costs and centers to the type of X, similar to wh…
tlienart Jan 28, 2019
192694a
Update src/kmeans.jl
alyst Jan 29, 2019
603422d
reverting the rounding of maxiter as per requested change
tlienart Jan 29, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
28 changes: 18 additions & 10 deletions src/kmeans.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -250,13 +250,15 @@ function update_centers!(
@inbounds for j = 1:n
cj = assignments[j]
1 <= cj <= k || error("assignment out of boundary.")

if to_update[cj]
xj = view(X, :, j)
if cweights[cj] > 0
centers[:, cj] .+= xj
for i = 1:d
tlienart marked this conversation as resolved.
Show resolved Hide resolved
centers[i, cj] += X[i, j]
end
else
centers[:, cj] .= xj
for i = 1:d
tlienart marked this conversation as resolved.
Show resolved Hide resolved
centers[i, cj] = X[i, j]
end
end
cweights[cj] += 1
end
Expand All @@ -265,7 +267,9 @@ function update_centers!(
# sum ==> mean
@inbounds for j = 1:k
if to_update[j]
centers[:, j] .*= 1 / cweights[j]
for i = 1:d
centers[i, j] /= cweights[j]
end
end
end
end
Expand Down Expand Up @@ -295,13 +299,15 @@ function update_centers!(
if wj > 0
cj = assignments[j]
1 <= cj <= k || error("assignment out of boundary.")

if to_update[cj]
xj = view(X, :, j)
if cweights[cj] > 0
centers[:, cj] .+= xj * wj
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess this * (plus the one below) was the reason for regression in weighted case as it creates the new vector.
.+= xj .* wj should avoid that.
I'm a little bit hesitant about converting everything to explicit loops as it defeats the purpose of having broadcast support in the language.

Copy link
Contributor Author

@tlienart tlienart Jan 28, 2019 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm happy with opening another PR for this kind of stuff (with tests & benchmarks) after this one gets merged so that we can isolate changes related to the initial problem from other changes?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

+1. Then I'd rather not introduce explicit loop conversions in this PR to keep it as small as possible.

Copy link
Contributor Author

@tlienart tlienart Jan 28, 2019 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

(FWIW, I just tried exactly what you said and it increases the allocation by a lot and doubles the time. There may be better things to do with copyto! etc but I must admit I'd rather get this PR done with first...)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought centers[:, cj] .+= calls setindex!() directly, but it actually creates centers[:, cj] vector behind the scenes (centers[:, cj] .= doesn't). So the better version is cej = view(centers, :, cj); cej .+= xj .* wj.

for i = 1:d
centers[i, cj] += X[i, j] * wj
end
else
centers[:, cj] .= xj * wj
for i = 1:d
centers[i, cj] = X[i, j] * wj
end
end
cweights[cj] += wj
end
Expand All @@ -311,7 +317,9 @@ function update_centers!(
# sum ==> mean
@inbounds for j = 1:k
if to_update[j]
centers[:, j] .*= 1 / cweights[j]
for i = 1:d
centers[i, j] /= cweights[j]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not 100% sure, but back in the days / was more expensive than *.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

with /= --> 132.276ms; 76.272ms
with *= 1 / ... --> 132.628ms; 77.268ms

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Note that the current code computes 1/cweights[j] only once, and then does the multiplication for each i. IIUC the timings for *= 1 / recompute the division it for each i, which isn't faster.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually, I think LLVM is smart enough to calculate 1/wj only once (at least that's what I see in @code_llvm for a toy example). But I still prefer the broadcast version (cej = view(centers, :, j); cej .*= 1/cweights[j]), where this would be definitely done once. :)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree, the broadcasted version was nicer. Better use view. Same above, where view was used but you removed it: use xj .* wj and everything should be OK.

end
end
end
end
Expand Down