Brachylog, 28 bytes

⟦₁⊇ᶠk⊇pSs₂ᶠ{c≠&⟨+/l⟩ᵐ<ᵈ}ᵐ∧Sl 

Try it online!

This is really damn slow. Takes about 30 seconds for N = 3, and it didn't complete after 12 minutes for N = 4.

Explanation

⟦₁ Take the range [1, …, Input] ⊇ᶠk Find all ordered subsets of that range, minus the empty set ⊇ Take an ordered subset of these subsets pS Take a permutation of that subset and call it S Ss₂ᶠ Find all substrings of 2 consecutive elements in S { }ᵐ Map for each of these substrings: c≠ All elements from both sets must be different &⟨+/l⟩ᵐ And the average of both sets (⟨xyz⟩ is a fork like in APL) <ᵈ Must be in strictly increasing order ∧Sl If all of this succeeds, the output is the length of L. 

Faster version, 39 bytes

⟦₁⊇ᶠk⊇{⟨+/l⟩/₁/₁}ᵒSs₂ᶠ{c≠&⟨+/l⟩ᵐ<₁}ᵐ∧Sl 

This takes about 50 seconds on my computer for N = 4.

This is the same program except we sort the subset of subsets by average instead of taking a random permutation. So we use {⟨+/l⟩/₁/₁}ᵒ instead of p.

{ }ᵒ Order by: ⟨+/l⟩ Average (fork of sum-divide-length) /₁/₁ Invert the average twice; this is used to get a float average 

We need to get a float average because I just discovered a ridiculous bug in which floats and integers don't compare by value but by type with ordering predicates (this is also why I use <ᵈ and not <₁ to compare both averages; the latter would require that double inversion trick to work).

CJam (81 bytes)

{YY@#(#{{2bW%ee{)*~}%}:Z~{)Z__1b\,d/\a+}%$}%{_,1>{2ew{z~~&!\~=>}%0&!}{,}?},:,:e>} 

Online demo. It should execute for input 4 in a reasonable time, but I wouldn't try it with higher inputs.

Dissection

{ e# Declare a block (anonymous function) YY@#(# e# There are 2^N subsets of [0, N), but the empty subset is problematic e# so we calculate 2^(2^N - 1) subsets of the non-empty subsets { e# Map integer to subset of non-empty subsets: { e# Define a block to map an bitset to its set indices; e.g. 9 => [0 3] 2bW%ee e# Convert to base 2, reverse, and index {)*~}% e# If the bit was set, keep the index }:Z e# Assign the block to variable Z ~ e# Evaluate it { e# Map those indices to non-empty subsets of [0, N): )Z e# Increment (to skip the empty set) and apply Z __1b\,d/ e# Sum one copy, take length of another, divide for average \a+ e# Wrap the subset and prepend its average value }% $ e# Sort (lexicographically, so by average value) }% { e# Filter out subsets of subsets with conflicts: _,1>{ e# If the length is greater than 1 2ew e# Take each consecutive pair of subsets { e# Map: z~ e# Zip and expand to get [score1 score2] [subset1 subset2] ~&!\ e# No element in common => 1 ~= e# Different scores => 0 > e# 1 iff both constraints are met }% 0&! e# 1 iff no consecutive pair failed the test }{ , e# Otherwise filter to those of length 1 }? }, :,:e> e# Map to size of subset and take the greatest } 

JavaScript (ES6), 175 bytes

A naive and rather slow recursive search. Takes about 15 seconds to compute the 7 first terms on TIO.

n=>(a=[...Array(n)].reduce(a=>[...a,...a.map(y=>[x,...y],x=n--)],[[]]),g=(p,b=a,n)=>a.map(a=>(m=a.map(n=>s+=++k*b.includes(n)?g:n,s=k=0)&&s/k)>p&&g(m,a,-~n),r=r>n?r:n))(r=0)|r 

Try it online!

or test this modified version that outputs the longest ex-increasing set sequence.

How?

We first compute the powerset of \$\{1,2,\dots,n\}\$ and store it in \$a\$:

a = [...Array(n)].reduce(a => [...a, ...a.map(y => [x, ...y], x = n--)], [[]] ) 

Recursive part:

g = ( // g = recursive function taking: p, // p = previous mean average b = a, // b = previous set n // n = sequence length ) => // a.map(a => // for each set a[] in a[]: (m = a.map(n => // for each value n in a[]: s += // update s: ++k * b.includes(n) ? // increment k; if n exists in b[]: g // invalidate the result (string / integer --> NaN) : // else: n, // add n to s s = k = 0) // start with s = k = 0; end of inner map() && s / k // m = s / k = new mean average ) > p // if it's greater than the previous one, && g(m, a, -~n), // do a recursive call with (m, a, n + 1) r = r > n ? r : n // keep track of the greatest length in r = max(r, n) ) // end of outer map() 

Python 3, 205 197 184 182 bytes

• Saved eight twenty-one bytes thanks to ovs.
• Saved two bytes thanks to ceilingcat.

f=lambda N,c=[[1]]:max([len(c)]+[f(N,c+[n])for k in range(N)for n in combinations(range(1,N+1),k+1)if not{*c[-1]}&{*n}and sum(c[-1])/len(c[-1])<sum(n)/len(n)]);from itertools import* 

Try it online!