# Efficient subsetting and averaging

The following code aims to average 4 groups of the values vector. Which elements belong to each group is determined by the indic vector and two split variables x1, x2:

• Group 1: those elements where indic is between 0 and x1,
• Group 2: those elements where indic is between x1 and x2,

and so on. The correspondence between the values vector and the indic is by position in the vector.

As a result, the code returns a list of the same length as values with the corresponding averages.

(* Input variables *)
indic = RandomChoice[Range[1,100], 1000];
values= RandomReal[{0,10}, 1000];

(* given certain breaks on indic, I want to calculate average on the values *)
myfun[x1_, x2_]:= Module[{avgVal, breaks, setAvg},
avgVal=ConstantArray[Missing, Length[values]];
breaks={0, x1,x2, Max[indic]};
setAvg[i_]:=Module[{idx, a},
idx=Position[indic, x_/;(x>breaks[[i]]&&x<=breaks[[i+1]])];
a=Mean[Extract[values,idx]];
avgVal=ReplacePart[avgVal, idx -> a]
];
setAvg /@ Range[1,3];
Return[avgVal];
];

myfun[20,50] (* averages values where indic < 20, between 20 and 50 and >50 *)


Minimal example:

indic  = {2,3,5,6,3,6,7,4};
values = Range[3,10];
x1=3;
x2=6;


Hence group 1 is the first, second, fifth element, group 2 comprises the third and the last element, and group 3 are the remaining elements.

The mean of group 1 is Mean[3,4,7]= 14/3, of group 2 is Mean[5,10]=15/2 and of group 3 the mean is Mean[6,8,9]=23/3. That is why myfun[x1,x2] should return

{14/3, 14/3, 15/2, 23/3, 14/3, 23/3, 23/3, 15/2}


I would like to discuss how to make this code faster. Is there a more Mathematica-like way to do it? Especially, I suspect the ReplacePart calls to be slow. Is there a way around it?

Any hint appreciated!

• Why don´t you describe in more detail what you want to do? Having to reverse-engineer your code is not so attractive. Commented Nov 27, 2014 at 12:13
• Sorry, there were errors in the code. I also edited the problem description. Commented Nov 27, 2014 at 12:34
• There may be a more efficient way, but is this on the right track?: {Pick[values, indic, x_ /; 0 < x < 3], Pick[values, indic, x_ /; 3 < x < 6]}. It's unclear what to do about the boundaries x1, x2, and what to do with the result, without the reverse-engineering Yves mentions. Commented Nov 27, 2014 at 12:37
• Thanks for mentioning Pick have not heard of it until now. The problem for me is to write the average vector. Is there a way to speed this up? Commented Nov 27, 2014 at 12:54

indic = RandomChoice[Range[1, 10], 100];
values = RandomReal[{0, 10}, 100];
xIntervals = {0, 3, 5, 8, Infinity};

pos = Position[indic, x_ /; (#1 <= x < #2)] & @@@ Partition[xIntervals, 2, 1];
l = Extract[values, #] & /@ pos;
means = Replace[l, {x__} :> ConstantArray[Mean@{x}, Length@{x}], {1}];
values1 = values;
(values1[[#1]] = #2) & @@@ Transpose[{Flatten@pos, Flatten@means}];


and now you have your desired result in values1

• You need to take care of the comparison limits (#1 <= x <= #2). I don't know how you want to handle the extremes Commented Nov 27, 2014 at 18:07
• This is about 15% faster than my code. I never saw Replace in use before. Thank you! Commented Nov 27, 2014 at 19:38
• @KarstenW. I didn't care about speed because you mentioned 1000 items. Perhaps I can improve it Commented Nov 27, 2014 at 19:43
• @KarstenW. And I'm stupid and haven't read your last sentence Commented Nov 27, 2014 at 19:45

I don't think I understand what you're trying to accomplish. Especially the way you use indic to set the outer bounds. For instance, note that it's not guaranteed that Max[indic] is larger than x2.

Anyway, the following might nonetheless be of use:

values = RandomReal[{0, 10}, 6];

myfun[n_, m_][values_] := Module[{split},
split = {values[[;; n]], values[[n + 1 ;; m]], values[[m + 1 ;;]]};
Flatten[ConstantArray[Mean[#], Length[#]] & /@ split]
];

myfun[2, 4][values]


{5.46557, 5.46557, 8.48547, 8.48547, 7.14661, 7.14661}

Note that your output will be different due to randomness of values.

Edit: If it's ok to have a balanced sampling, you could try random sampling as follows:

values = RandomReal[{0, 10}, 10];
range = Range[Length[values]];

(* get three sample with balanced sample size *)
numSamples = 3;
size = Ceiling[Length[range]/numSamples];
samples = Partition[RandomSample[range], size, size, {1, 1}, {}];

(* define a dictionary that associates each index with an appropriate mean *)
means =
Rule @@ # & /@
Transpose@{#, ConstantArray[Mean@values[[#]], Length@values[[#]]]} & /@
samples // Association;

(* output *)
means /@ range


{5.24007, 4.63392, 5.24007, 5.5112, 4.63392, 5.5112}

Let me know if I need to unpack the definition for means.

• I'd like a vector of the same length as values, how can I do this? Commented Nov 27, 2014 at 13:04
• @KarstenW. I edited the answer. Is this closer to what you were looking for?
– Kris
Commented Nov 27, 2014 at 14:29
• The output looks now as desired, however I do not really understand the samples variable. Commented Nov 27, 2014 at 19:41
• samples is a partitioned version of the randomly shuffled version of range, where range contains the indices for values. In other words, Sort[Join @@ samples] gives back range.
– Kris
Commented Nov 27, 2014 at 21:58