diff --git a/base/array.jl b/base/array.jl
index fd38e93d19b38..eb7bba7af66fc 100644
--- a/base/array.jl
+++ b/base/array.jl
@@ -9,6 +9,14 @@ typealias StridedVector{T,A<:Array}   Union(Vector{T} , SubArray{T,1,A})
 typealias StridedMatrix{T,A<:Array}   Union(Matrix{T} , SubArray{T,2,A})
 typealias StridedVecOrMat{T} Union(StridedVector{T}, StridedMatrix{T})
 
+type Dim{A<:AbstractVector{Int}}
+    region::A
+end
+Dim(i::Int) = Dim([i])
+Dim(A::AbstractArray{Int}) = Dim(vec(A))
+Dim(R::Range) = Dim(vec(R))
+Dim(A...) = Dim([A...])
+
 ## Basic functions ##
 
 size(a::Array) = arraysize(a)
@@ -1425,19 +1433,19 @@ function reducedim(f::Function, A, region, v0, R)
 end
 end
 
-max{T}(A::AbstractArray{T}, b::(), region) = reducedim(max,A,region,typemin(T))
-min{T}(A::AbstractArray{T}, b::(), region) = reducedim(min,A,region,typemax(T))
-sum{T}(A::AbstractArray{T}, region)  = reducedim(+,A,region,zero(T))
-prod{T}(A::AbstractArray{T}, region) = reducedim(*,A,region,one(T))
+max{T}(A::AbstractArray{T}, d::Dim) = reducedim(max,A,d.region,typemin(T))
+min{T}(A::AbstractArray{T}, d::Dim) = reducedim(min,A,d.region,typemax(T))
+sum{T}(A::AbstractArray{T}, d::Dim)  = reducedim(+,A,d.region,zero(T))
+prod{T}(A::AbstractArray{T}, d::Dim) = reducedim(*,A,d.region,one(T))
 
-all(A::AbstractArray{Bool}, region) = reducedim(all,A,region,true)
-any(A::AbstractArray{Bool}, region) = reducedim(any,A,region,false)
-sum(A::AbstractArray{Bool}, region) = reducedim(+,A,region,0,similar(A,Int,reduced_dims(A,region)))
-sum(A::AbstractArray{Bool}) = sum(A, [1:ndims(A)])[1]
-sum(A::StridedArray{Bool})  = sum(A, [1:ndims(A)])[1]
+all(A::AbstractArray{Bool}, d::Dim) = reducedim(all,A,d.region,true)
+any(A::AbstractArray{Bool}, d::Dim) = reducedim(any,A,d.region,false)
+sum(A::AbstractArray{Bool}, d::Dim) = reducedim(+,A,d.region,0,similar(A,Int,reduced_dims(A,d.region)))
+sum(A::AbstractArray{Bool}) = sum(A, Dim(1:ndims(A)) )[1]
+sum(A::StridedArray{Bool})  = sum(A, Dim(1:ndims(A)) )[1]
 prod(A::AbstractArray{Bool}) =
     error("use all() instead of prod() for boolean arrays")
-prod(A::AbstractArray{Bool}, region) =
+prod(A::AbstractArray{Bool}, d::Dim) =
     error("use all() instead of prod() for boolean arrays")
 
 function sum{T}(A::StridedArray{T})
diff --git a/base/bitarray.jl b/base/bitarray.jl
index d0d33b17b921d..5eeed9443a327 100644
--- a/base/bitarray.jl
+++ b/base/bitarray.jl
@@ -1604,7 +1604,7 @@ nonzeros(B::BitArray) = trues(nnz(B))
 
 ## Reductions ##
 
-sum(A::BitArray, region) = reducedim(+,A,region,0,Array(Int,reduced_dims(A,region)))
+sum(A::BitArray, d::Dim) = reducedim(+,A,d.region,0,Array(Int,reduced_dims(A,d.region)))
 
 sum(B::BitArray) = nnz(B)
 
diff --git a/base/exports.jl b/base/exports.jl
index 3c32c971758ac..972541d2e0b85 100644
--- a/base/exports.jl
+++ b/base/exports.jl
@@ -38,6 +38,7 @@ export
     ComplexPair,
     DArray,
     Dict,
+    Dim,
     Dims,
     EachLine,
     Enumerate,
diff --git a/base/linalg.jl b/base/linalg.jl
index e742a5be32ee8..51a57e7d8bf4b 100644
--- a/base/linalg.jl
+++ b/base/linalg.jl
@@ -60,11 +60,11 @@ function norm(A::AbstractMatrix, p::Number)
     elseif m == 1 || n == 1
         a = norm(reshape(A, length(A)), p)
     elseif p == 1
-        a = max(sum(abs(A),1))
+        a = max(sum(abs(A),Dim(1)))
     elseif p == 2
         a = max(svdvals(A))
     elseif p == Inf
-        a = max(sum(abs(A),2))
+        a = max(sum(abs(A),Dim(2)))
     else
         error("invalid parameter p given to compute matrix norm")
     end
diff --git a/base/sparse.jl b/base/sparse.jl
index 42a42a09fca36..def8f4cf9621b 100644
--- a/base/sparse.jl
+++ b/base/sparse.jl
@@ -680,7 +680,7 @@ end # macro
 
 # TODO: Should the results of sparse reductions be sparse?
 function reducedim{Tv,Ti}(f::Function, A::SparseMatrixCSC{Tv,Ti}, region, v0)
-    if region == 1
+    if region == [1]
 
         S = Array(Tv, 1, A.n)
         for i = 1 : A.n
@@ -695,7 +695,7 @@ function reducedim{Tv,Ti}(f::Function, A::SparseMatrixCSC{Tv,Ti}, region, v0)
         end
         return S
 
-    elseif region == 2
+    elseif region == [2]
 
         S = fill(v0, A.m, 1)
         rcounts = zeros(Ti, A.m)
@@ -709,7 +709,7 @@ function reducedim{Tv,Ti}(f::Function, A::SparseMatrixCSC{Tv,Ti}, region, v0)
         end
         return S
 
-    elseif region == (1,2)
+    elseif region == [1,2]
 
         S = v0
         for i = 1 : A.n, j = A.colptr[i] : A.colptr[i+1]-1
@@ -726,17 +726,17 @@ function reducedim{Tv,Ti}(f::Function, A::SparseMatrixCSC{Tv,Ti}, region, v0)
     end
 end
 
-max{T}(A::SparseMatrixCSC{T}) = reducedim(max,A,(1,2),typemin(T))
-max{T}(A::SparseMatrixCSC{T}, b::(), region) = reducedim(max,A,region,typemin(T))
+max{T}(A::SparseMatrixCSC{T}) = reducedim(max,A,[1,2],typemin(T))
+max{T}(A::SparseMatrixCSC{T}, d::Dim) = reducedim(max,A,d.region,typemin(T))
 
-min{T}(A::SparseMatrixCSC{T}) = reducedim(min,A,(1,2),typemax(T))
-min{T}(A::SparseMatrixCSC{T}, b::(), region) = reducedim(min,A,region,typemax(T))
+min{T}(A::SparseMatrixCSC{T}) = reducedim(min,A,[1,2],typemax(T))
+min{T}(A::SparseMatrixCSC{T}, d::Dim) = reducedim(min,A,d.region,typemax(T))
 
-sum{T}(A::SparseMatrixCSC{T}) = reducedim(+,A,(1,2),zero(T))
-sum{T}(A::SparseMatrixCSC{T}, region)  = reducedim(+,A,region,zero(T))
+sum{T}(A::SparseMatrixCSC{T}) = reducedim(+,A,[1,2],zero(T))
+sum{T}(A::SparseMatrixCSC{T}, d::Dim)  = reducedim(+,A,d.region,zero(T))
 
-prod{T}(A::SparseMatrixCSC{T}) = reducedim(*,A,(1,2),one(T))
-prod{T}(A::SparseMatrixCSC{T}, region) = reducedim(*,A,region,one(T))
+prod{T}(A::SparseMatrixCSC{T}) = reducedim(*,A,[1,2],one(T))
+prod{T}(A::SparseMatrixCSC{T}, d::Dim) = reducedim(*,A,d.region,one(T))
 
 #all(A::SparseMatrixCSC{Bool}, region) = reducedim(all,A,region,true)
 #any(A::SparseMatrixCSC{Bool}, region) = reducedim(any,A,region,false)
diff --git a/base/statistics.jl b/base/statistics.jl
index 31674c68d3795..3c8f775a3d6f1 100644
--- a/base/statistics.jl
+++ b/base/statistics.jl
@@ -12,6 +12,7 @@ function mean(iterable)
     end
     return total/count
 end
+mean(v::AbstractArray, d::Dim) = sum(v, d) / prod(size(v)[d.region])
 
 function median!{T<:Real}(v::AbstractVector{T})
     isempty(v) && error("median of an empty array is undefined")
@@ -20,6 +21,7 @@ function median!{T<:Real}(v::AbstractVector{T})
     isodd(length(v)) ? float(v[div(end+1,2)]) : (v[div(end,2)]+v[div(end,2)+1])/2
 end
 median{T<:Real}(v::AbstractArray{T}) = median!(copy(reshape(v, length(v))))
+#median{T<:Real}(v::AbstractArray{T}, d::Dim) = # ?? TODO
 
 ## variance with known mean
 function var(v::AbstractVector, m::Number, corrected::Bool)
@@ -33,6 +35,7 @@ end
 var(v::AbstractVector, m::Number) = var(v, m, true)
 var(v::AbstractArray, m::Number, corrected::Bool) = var(reshape(v, length(v)), m, corrected)
 var(v::AbstractArray, m::Number) = var(v, m, true)
+#var(v::AbstractArray, m::Number, d::Dim) = # ?? TODO
 function var(v::Ranges, m::Number, corrected::Bool)
     f = first(v) - m
     s = step(v)
@@ -59,15 +62,20 @@ function var(v::Ranges, corrected::Bool)
 end
 var(v::AbstractVector, corrected::Bool) = var(v, mean(v), corrected)
 var(v::AbstractArray, corrected::Bool) = var(reshape(v, length(v)), corrected)
+#var(v::AbstractArray, d::Dim, corrected::Bool) = # ?? TODO
+#var(v::AbstractArray, d::Dim) = var(v, d, true)
 var(v::AbstractArray) = var(v, true)
 
 ## standard deviation with known mean
 std(v::AbstractArray, m::Number, corrected::Bool) = sqrt(var(v, m, corrected))
 std(v::AbstractArray, m::Number) = std(v, m, true)
+# std(v::AbstractArray, m::Number, d::Dim) = sqrt(var(v, m, d, true))
+# std(v::AbstractArray, m::Number, d::Dim, corrected::Bool) = sqrt(var(v, m, d, corrected))
 
 ## standard deviation
 std(v::AbstractArray, corrected::Bool) = std(v, mean(v), corrected)
 std(v::AbstractArray) = std(v, true)
+# std(v::AbstractArray, d::Dim) = std(v, mean(v, d), d, corrected)
 std(v::Ranges, corrected::Bool) = sqrt(var(v, corrected))
 std(v::Ranges) = std(v, true)
 
diff --git a/test/arrayops.jl b/test/arrayops.jl
index 2f9863792bb02..25b15906aab42 100644
--- a/test/arrayops.jl
+++ b/test/arrayops.jl
@@ -219,7 +219,7 @@ for i=1:16
     z[i] = i
 end
 
-@test sum(z) == sum(z,(1,2,3,4))[1] == 136
+@test sum(z) == sum(z,Dim(1,2,3,4))[1] == 136
 
 v = cell(2,2,1,1)
 v[1,1,1,1] = 28.0
@@ -227,7 +227,7 @@ v[1,2,1,1] = 36.0
 v[2,1,1,1] = 32.0
 v[2,2,1,1] = 40.0
 
-@test isequal(v,sum(z,(3,4)))
+@test isequal(v,sum(z,Dim(3,4)))
 
 ## large matrices transpose ##
 
diff --git a/test/sparse.jl b/test/sparse.jl
index 0622e623e87c0..65f9bc89e79d4 100644
--- a/test/sparse.jl
+++ b/test/sparse.jl
@@ -63,8 +63,8 @@ end
 
 # reductions
 @test sum(se33)[1] == 3.0
-@test sum(se33, 1) == [1.0 1.0 1.0]
-@test sum(se33, 2) == [1.0 1.0 1.0]'
+@test sum(se33, Dim(1)) == [1.0 1.0 1.0]
+@test sum(se33, Dim(2)) == [1.0 1.0 1.0]'
 @test prod(se33)[1] == 0.0
-@test prod(se33, 1) == [0.0 0.0 0.0]
-@test prod(se33, 2) == [0.0 0.0 0.0]'
+@test prod(se33, Dim(1)) == [0.0 0.0 0.0]
+@test prod(se33, Dim(2)) == [0.0 0.0 0.0]'